/[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

revision 472 by jgs, Fri Jan 27 01:50:59 2006 UTC revision 967 by gross, Tue Feb 13 09:40:12 2007 UTC
# Line 13  Line 13 
13   ******************************************************************************   ******************************************************************************
14  */  */
15    
16    #ifdef PASO_MPI
17    #include <mpi.h>
18    #endif
19  #include "MeshAdapter.h"  #include "MeshAdapter.h"
20    
21    #include "escript/Data.h"
22    #include "escript/DataFactory.h"
23    
24  using namespace std;  using namespace std;
25  using namespace escript;  using namespace escript;
26    
# Line 67  Finley_Mesh* MeshAdapter::getFinley_Mesh Line 73  Finley_Mesh* MeshAdapter::getFinley_Mesh
73    
74  void MeshAdapter::write(const std::string& fileName) const  void MeshAdapter::write(const std::string& fileName) const
75  {  {
76    char fName[fileName.size()+1];    char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
77    strcpy(fName,fileName.c_str());    strcpy(fName,fileName.c_str());
78    Finley_Mesh_write(m_finleyMesh.get(),fName);    Finley_Mesh_write(m_finleyMesh.get(),fName);
79    checkFinleyError();    checkFinleyError();
80      TMPMEMFREE(fName);
81  }  }
82    
 // void MeshAdapter::getTagList(int functionSpaceType,  
 //                  int* numTags) const  
 // {  
 //   Finley_Mesh_tagList(m_finleyMesh.get(),functionSpaceType,tagList,numTags);  
 //   return;  
 // }  
   
83  string MeshAdapter::getDescription() const  string MeshAdapter::getDescription() const
84  {  {
85    return "FinleyMesh";    return "FinleyMesh";
# Line 164  int MeshAdapter::getDiracDeltaFunctionCo Line 164  int MeshAdapter::getDiracDeltaFunctionCo
164  }  }
165    
166  //  //
 // 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;  
   }  
   if (*tagList==NULL) {  
     stringstream temp;  
     temp << "Error - no tags available for " << functionSpaceType << " for domain: " << getDescription();  
     throw FinleyAdapterException(temp.str());  
   }  
   return;  
 }  
   
 //  
 // returns a pointer to the reference no list of samples of functionSpaceType  
 //  
 void MeshAdapter::getReferenceNoList(int functionSpaceType, int** referenceNoList,  
                  int* numReferenceNo) const  
 {  
   *referenceNoList=NULL;  
   *numReferenceNo=0;  
   Finley_Mesh* mesh=m_finleyMesh.get();  
   switch (functionSpaceType) {  
   case(Nodes):  
     if (mesh->Nodes!=NULL) {  
       *referenceNoList=mesh->Nodes->Id;  
       *numReferenceNo=mesh->Nodes->numNodes;  
     }  
     break;  
   case(Elements):  
     if (mesh->Elements!=NULL) {  
       *referenceNoList=mesh->Elements->Id;  
       *numReferenceNo=mesh->Elements->numElements;  
     }  
     break;  
   case(FaceElements):  
     if (mesh->FaceElements!=NULL) {  
       *referenceNoList=mesh->FaceElements->Id;  
       *numReferenceNo=mesh->FaceElements->numElements;  
     }  
     break;  
   case(Points):  
     if (mesh->Points!=NULL) {  
       *referenceNoList=mesh->Points->Id;  
       *numReferenceNo=mesh->Points->numElements;  
     }  
     break;  
   case(ContactElementsZero):  
     if (mesh->ContactElements!=NULL) {  
       *referenceNoList=mesh->ContactElements->Id;  
       *numReferenceNo=mesh->ContactElements->numElements;  
     }  
     break;  
   case(ContactElementsOne):  
     if (mesh->ContactElements!=NULL) {  
       *referenceNoList=mesh->ContactElements->Id;  
       *numReferenceNo=mesh->ContactElements->numElements;  
     }  
     break;  
   case(DegreesOfFreedom):  
     if (mesh->Nodes!=NULL) {  
       *referenceNoList=NULL;  
       *numReferenceNo=0;  
     }  
     break;  
   case(ReducedDegreesOfFreedom):  
     if (mesh->Nodes!=NULL) {  
       *referenceNoList=NULL;  
       *numReferenceNo=0;  
     }  
     break;  
   default:  
     stringstream temp;  
     temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();  
     throw FinleyAdapterException(temp.str());  
     break;  
   }  
   if (*referenceNoList==NULL) {  
     stringstream temp;  
     temp << "Error - reference number list available for " << functionSpaceType << " for domain: " << getDescription();  
     throw FinleyAdapterException(temp.str());  
   }  
   return;  
 }  
   
 //  
167  // return the spatial dimension of the Mesh:  // return the spatial dimension of the Mesh:
168  //  //
169  int MeshAdapter::getDim() const  int MeshAdapter::getDim() const
# Line 364  pair<int,int> MeshAdapter::getDataShape( Line 220  pair<int,int> MeshAdapter::getDataShape(
220        case(DegreesOfFreedom):        case(DegreesOfFreedom):
221             if (mesh->Nodes!=NULL) {             if (mesh->Nodes!=NULL) {
222               numDataPointsPerSample=1;               numDataPointsPerSample=1;
223    #ifndef PASO_MPI
224               numSamples=mesh->Nodes->numDegreesOfFreedom;               numSamples=mesh->Nodes->numDegreesOfFreedom;
225    #else
226                 numSamples=mesh->Nodes->degreeOfFreedomDistribution->numLocal;
227    #endif
228             }             }
229             break;             break;
230        case(ReducedDegreesOfFreedom):        case(ReducedDegreesOfFreedom):
231             if (mesh->Nodes!=NULL) {             if (mesh->Nodes!=NULL) {
232               numDataPointsPerSample=1;               numDataPointsPerSample=1;
233    #ifndef PASO_MPI
234               numSamples=mesh->Nodes->reducedNumDegreesOfFreedom;               numSamples=mesh->Nodes->reducedNumDegreesOfFreedom;
235    #else
236                 numSamples=mesh->Nodes->reducedDegreeOfFreedomDistribution->numLocal;
237    #endif
238             }             }
239             break;             break;
240        default:        default:
# Line 386  pair<int,int> MeshAdapter::getDataShape( Line 250  pair<int,int> MeshAdapter::getDataShape(
250  // 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:
251  //  //
252  void MeshAdapter::addPDEToSystem(  void MeshAdapter::addPDEToSystem(
253                       SystemMatrixAdapter& mat, Data& rhs,                       SystemMatrixAdapter& mat, escript::Data& rhs,
254                       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,
255                       const Data& d, const Data& y,                       const escript::Data& d, const escript::Data& y,
256                       const Data& d_contact,const Data& y_contact) const                       const escript::Data& d_contact,const escript::Data& y_contact) const
257  {  {
258       escriptDataC _rhs=rhs.getDataC();
259       escriptDataC _A  =A.getDataC();
260       escriptDataC _B=B.getDataC();
261       escriptDataC _C=C.getDataC();
262       escriptDataC _D=D.getDataC();
263       escriptDataC _X=X.getDataC();
264       escriptDataC _Y=Y.getDataC();
265       escriptDataC _d=d.getDataC();
266       escriptDataC _y=y.getDataC();
267       escriptDataC _d_contact=d_contact.getDataC();
268       escriptDataC _y_contact=y_contact.getDataC();
269    
270     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
271     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(),&(rhs.getDataC()),  
272                         &(A.getDataC()),&(B.getDataC()),&(C.getDataC()),&(D.getDataC()),&(X.getDataC()),&(Y.getDataC()));     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(), &_rhs, &_A, &_B, &_C, &_D, &_X, &_Y );
273     checkFinleyError();     checkFinleyError();
274     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,  
275                    mat.getPaso_SystemMatrix(),     Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, mat.getPaso_SystemMatrix(), &_rhs, 0, 0, 0, &_d, 0, &_y );
                   &(rhs.getDataC()),  
                                   &(d.getDataC()),&(y.getDataC()),  
                                   Finley_Assemble_handelShapeMissMatch_Mean_out);  
276     checkFinleyError();     checkFinleyError();
277     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->ContactElements,  
278                    mat.getPaso_SystemMatrix(),     Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, mat.getPaso_SystemMatrix(), &_rhs , 0, 0, 0, &_d_contact, 0, &_y_contact );
                   &(rhs.getDataC()),  
                                   &(d_contact.getDataC()),  
                   &(y_contact.getDataC()),  
                                   Finley_Assemble_handelShapeMissMatch_Step_out);  
279     checkFinleyError();     checkFinleyError();
280  }  }
281    
282  //  //
283  // adds linear PDE of second order into the right hand side only  // adds linear PDE of second order into the right hand side only
284  //  //
285  void MeshAdapter::addPDEToRHS( Data& rhs,  void MeshAdapter::addPDEToRHS( escript::Data& rhs, const  escript::Data& X,const  escript::Data& Y, const escript::Data& y, const escript::Data& y_contact) const
                      const  Data& X,const  Data& Y, const Data& y, const Data& y_contact) const  
286  {  {
287     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
288    
289     // Finley_Assemble_PDE_RHS(mesh->Nodes,mesh->Elements,&(rhs.getDataC()),&(X.getDataC()),&(Y.getDataC()));     escriptDataC _rhs=rhs.getDataC();
290     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,0,&(rhs.getDataC()),0,0,0,0,&(X.getDataC()),&(Y.getDataC()));     escriptDataC _X=X.getDataC();
291     checkFinleyError();     escriptDataC _Y=Y.getDataC();
292       escriptDataC _y=y.getDataC();
293       escriptDataC _y_contact=y_contact.getDataC();
294    
295     // Finley_Assemble_RobinCondition_RHS(mesh->Nodes,mesh->FaceElements,&(rhs.getDataC()),&(y.getDataC()),Finley_Assemble_handelShapeMissMatch_Mean_out);     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements, 0, &_rhs, 0, 0, 0, 0, &_X, &_Y );
296     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,0,&(rhs.getDataC()),0,&(y.getDataC()),Finley_Assemble_handelShapeMissMatch_Mean_out);     checkFinleyError();
297    
298       Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, 0, &_rhs, 0, 0, 0, 0, 0, &_y );
299     checkFinleyError();     checkFinleyError();
300     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->ContactElements,0,&(rhs.getDataC()),0,&(y_contact.getDataC()),Finley_Assemble_handelShapeMissMatch_Step_out);  
301     // Finley_Assemble_RobinCondition_RHS(mesh->Nodes,mesh->ContactElements,&(rhs.getDataC()),&(y_contact.getDataC()),Finley_Assemble_handelShapeMissMatch_Step_out);     Finley_Assemble_PDE(mesh->Nodes,mesh->ContactElements, 0, &_rhs , 0, 0, 0, 0, 0, &_y_contact );
302     checkFinleyError();     checkFinleyError();
303  }  }
304    
305  //  //
306  // interpolates data between different function spaces:  // interpolates data between different function spaces:
307  //  //
308  void MeshAdapter::interpolateOnDomain(Data& target,const Data& in) const  void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
309  {  {
310    const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());    const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());
311    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
312    if (inDomain!=*this)    if (inDomain!=*this)  
313       throw FinleyAdapterException("Error - Illegal domain of interpolant.");      throw FinleyAdapterException("Error - Illegal domain of interpolant.");
314    if (targetDomain!=*this)    if (targetDomain!=*this)
315       throw FinleyAdapterException("Error - Illegal domain of interpolation target.");      throw FinleyAdapterException("Error - Illegal domain of interpolation target.");
316    
317    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
318    switch(in.getFunctionSpace().getTypeCode()) {    switch(in.getFunctionSpace().getTypeCode()) {
# Line 505  void MeshAdapter::interpolateOnDomain(Da Line 376  void MeshAdapter::interpolateOnDomain(Da
376             break;             break;
377          }          }
378          break;          break;
379       case(DegreesOfFreedom):       case(DegreesOfFreedom):      
380          switch(target.getFunctionSpace().getTypeCode()) {          switch(target.getFunctionSpace().getTypeCode()) {
381             case(ReducedDegreesOfFreedom):             case(ReducedDegreesOfFreedom):
382             case(DegreesOfFreedom):             case(DegreesOfFreedom):
383             case(Nodes):             case(Nodes):
384                Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));                Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));
385                break;                break;
386    #ifndef PASO_MPI
387             case(Elements):             case(Elements):
388                Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(in.getDataC()),&(target.getDataC()));                Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(in.getDataC()),&(target.getDataC()));
389                break;                break;
# Line 525  void MeshAdapter::interpolateOnDomain(Da Line 397  void MeshAdapter::interpolateOnDomain(Da
397             case(ContactElementsOne):             case(ContactElementsOne):
398                Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));                Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
399               break;               break;
400    #else
401               /* need to copy Degrees of freedom data to nodal data so that the external values are available */
402               case(Elements):
403               {
404                  escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
405                  Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&(nodeTemp.getDataC()),&(target.getDataC()));
406                  break;
407               }
408               case(FaceElements):
409               {
410                  escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
411                  Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&(nodeTemp.getDataC()),&(target.getDataC()));
412                  break;
413               }
414               case(Points):
415               {
416                  escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
417                  Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&(nodeTemp.getDataC()),&(target.getDataC()));
418                  break;
419               }
420               case(ContactElementsZero):
421               case(ContactElementsOne):
422               {
423                  escript::Data nodeTemp( in, continuousFunction(asAbstractContinuousDomain()) );
424                  Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(nodeTemp.getDataC()),&(target.getDataC()));
425                 break;
426               }
427    #endif
428             default:             default:
429               stringstream temp;               stringstream temp;
430               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
# Line 575  void MeshAdapter::interpolateOnDomain(Da Line 475  void MeshAdapter::interpolateOnDomain(Da
475  //  //
476  // copies the locations of sample points into x:  // copies the locations of sample points into x:
477  //  //
478  void MeshAdapter::setToX(Data& arg) const  void MeshAdapter::setToX(escript::Data& arg) const
479  {  {
480    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
481    if (argDomain!=*this)    if (argDomain!=*this)
# Line 585  void MeshAdapter::setToX(Data& arg) cons Line 485  void MeshAdapter::setToX(Data& arg) cons
485    if (arg.getFunctionSpace().getTypeCode()==Nodes) {    if (arg.getFunctionSpace().getTypeCode()==Nodes) {
486       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(arg.getDataC()));       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(arg.getDataC()));
487    } else {    } else {
488       Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);       escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
489       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(tmp_data.getDataC()));       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(tmp_data.getDataC()));
490       // this is then interpolated onto arg:       // this is then interpolated onto arg:
491       interpolateOnDomain(arg,tmp_data);       interpolateOnDomain(arg,tmp_data);
# Line 596  void MeshAdapter::setToX(Data& arg) cons Line 496  void MeshAdapter::setToX(Data& arg) cons
496  //  //
497  // 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:
498  //  //
499  void MeshAdapter::setToNormal(Data& normal) const  void MeshAdapter::setToNormal(escript::Data& normal) const
500  {  {
501    const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());    const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());
502    if (normalDomain!=*this)    if (normalDomain!=*this)
# Line 637  void MeshAdapter::setToNormal(Data& norm Line 537  void MeshAdapter::setToNormal(Data& norm
537  //  //
538  // interpolates data to other domain:  // interpolates data to other domain:
539  //  //
540  void MeshAdapter::interpolateACross(Data& target,const Data& source) const  void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
541  {  {
542    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
543    if (targetDomain!=*this)    if (targetDomain!=*this)
# Line 649  void MeshAdapter::interpolateACross(Data Line 549  void MeshAdapter::interpolateACross(Data
549  //  //
550  // calculates the integral of a function defined of arg:  // calculates the integral of a function defined of arg:
551  //  //
552  void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const Data& arg) const  void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
553  {  {
554    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
555    if (argDomain!=*this)    if (argDomain!=*this)
# Line 693  void MeshAdapter::setToIntegrals(std::ve Line 593  void MeshAdapter::setToIntegrals(std::ve
593  //  //
594  // calculates the gradient of arg:  // calculates the gradient of arg:
595  //  //
596  void MeshAdapter::setToGradient(Data& grad,const Data& arg) const  void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
597  {  {
598    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
599    if (argDomain!=*this)    if (argDomain!=*this)
# Line 703  void MeshAdapter::setToGradient(Data& gr Line 603  void MeshAdapter::setToGradient(Data& gr
603       throw FinleyAdapterException("Error - Illegal domain of gradient");       throw FinleyAdapterException("Error - Illegal domain of gradient");
604    
605    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
606      escriptDataC nodeDataC;
607    #ifdef PASO_MPI
608      escript::Data nodeTemp( arg, continuousFunction(asAbstractContinuousDomain()) );
609      if( arg.getFunctionSpace().getTypeCode() != Nodes )
610      {
611        Finley_Assemble_CopyNodalData(mesh->Nodes,&(nodeTemp.getDataC()),&(arg.getDataC()));
612        nodeDataC = nodeTemp.getDataC();
613      }
614      else
615        nodeDataC = arg.getDataC();
616    #else
617      nodeDataC = arg.getDataC();
618    #endif
619    switch(grad.getFunctionSpace().getTypeCode()) {    switch(grad.getFunctionSpace().getTypeCode()) {
620         case(Nodes):         case(Nodes):
621            throw FinleyAdapterException("Error - Gradient at nodes is not supported.");            throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
622            break;            break;
623         case(Elements):         case(Elements):
624            Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&nodeDataC);
625            break;            break;
626         case(FaceElements):         case(FaceElements):
627            Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&nodeDataC);
628            break;            break;
629         case(Points):         case(Points):
630            throw FinleyAdapterException("Error - Gradient at points is not supported.");            throw FinleyAdapterException("Error - Gradient at points is not supported.");
631            break;            break;
632         case(ContactElementsZero):         case(ContactElementsZero):
633            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
634            break;            break;
635         case(ContactElementsOne):         case(ContactElementsOne):
636            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&nodeDataC);
637            break;            break;
638         case(DegreesOfFreedom):         case(DegreesOfFreedom):
639            throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");            throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
# Line 740  void MeshAdapter::setToGradient(Data& gr Line 653  void MeshAdapter::setToGradient(Data& gr
653  //  //
654  // returns the size of elements:  // returns the size of elements:
655  //  //
656  void MeshAdapter::setToSize(Data& size) const  void MeshAdapter::setToSize(escript::Data& size) const
657  {  {
658    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
659    escriptDataC tmp=size.getDataC();    escriptDataC tmp=size.getDataC();
# Line 777  void MeshAdapter::setToSize(Data& size) Line 690  void MeshAdapter::setToSize(Data& size)
690  }  }
691    
692  // sets the location of nodes:  // sets the location of nodes:
693  void MeshAdapter::setNewX(const Data& new_x)  void MeshAdapter::setNewX(const escript::Data& new_x)
694  {  {
695    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
696      escriptDataC tmp;
697    const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());    const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());
698    if (newDomain!=*this)    if (newDomain!=*this)
699       throw FinleyAdapterException("Error - Illegal domain of new point locations");       throw FinleyAdapterException("Error - Illegal domain of new point locations");
700    Finley_NodeFile_setCoordinates(mesh->Nodes,&(new_x.getDataC()));    tmp = new_x.getDataC();
701      Finley_Mesh_setCoordinates(mesh,&tmp);
702    checkFinleyError();    checkFinleyError();
703  }  }
704    
# Line 792  void MeshAdapter::saveDX(const std::stri Line 707  void MeshAdapter::saveDX(const std::stri
707  {  {
708      int MAX_namelength=256;      int MAX_namelength=256;
709      const int num_data=boost::python::extract<int>(arg.attr("__len__")());      const int num_data=boost::python::extract<int>(arg.attr("__len__")());
710      char names[num_data][MAX_namelength];    /* win32 refactor */
711      char* c_names[num_data];    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
712      escriptDataC data[num_data];    for(int i=0;i<num_data;i++)
713      escriptDataC* ptr_data[num_data];    {
714        names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
715      }
716    
717      char* *c_names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
718      escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
719      escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
720    
721      boost::python::list keys=arg.keys();      boost::python::list keys=arg.keys();
722      for (int i=0;i<num_data;++i) {      for (int i=0;i<num_data;++i) {
723           Data& d=boost::python::extract<Data&>(arg[keys[i]]);           escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
724           if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)           if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
725               throw FinleyAdapterException("Error  in saveVTK: Data must be defined on same Domain");               throw FinleyAdapterException("Error  in saveDX: Data must be defined on same Domain");
726           data[i]=d.getDataC();           data[i]=d.getDataC();
727           ptr_data[i]=&(data[i]);           ptr_data[i]=&(data[i]);
728           std::string n=boost::python::extract<std::string>(keys[i]);           std::string n=boost::python::extract<std::string>(keys[i]);
# Line 815  void MeshAdapter::saveDX(const std::stri Line 736  void MeshAdapter::saveDX(const std::stri
736      }      }
737      Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);      Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
738      checkFinleyError();      checkFinleyError();
739        
740          /* win32 refactor */
741      TMPMEMFREE(c_names);
742      TMPMEMFREE(data);
743      TMPMEMFREE(ptr_data);
744      for(int i=0;i<num_data;i++)
745      {
746        TMPMEMFREE(names[i]);
747      }
748      TMPMEMFREE(names);
749    
750      return;      return;
751  }  }
752    
# Line 823  void MeshAdapter::saveVTK(const std::str Line 755  void MeshAdapter::saveVTK(const std::str
755  {  {
756      int MAX_namelength=256;      int MAX_namelength=256;
757      const int num_data=boost::python::extract<int>(arg.attr("__len__")());      const int num_data=boost::python::extract<int>(arg.attr("__len__")());
758      char names[num_data][MAX_namelength];    /* win32 refactor */
759      char* c_names[num_data];    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
760      escriptDataC data[num_data];    for(int i=0;i<num_data;i++)
761      escriptDataC* ptr_data[num_data];    {
762        names[i] = (MAX_namelength>0) ? TMPMEMALLOC(MAX_namelength,char) : (char*)NULL;
763      }
764    
765      char* *c_names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
766      escriptDataC *data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC) : (escriptDataC*)NULL;
767      escriptDataC* *ptr_data = (num_data>0) ? TMPMEMALLOC(num_data,escriptDataC*) : (escriptDataC**)NULL;
768    
769      boost::python::list keys=arg.keys();      boost::python::list keys=arg.keys();
770      for (int i=0;i<num_data;++i) {      for (int i=0;i<num_data;++i) {
771           Data& d=boost::python::extract<Data&>(arg[keys[i]]);           escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
772           if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)           if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
773               throw FinleyAdapterException("Error  in saveVTK: Data must be defined on same Domain");               throw FinleyAdapterException("Error  in saveVTK: Data must be defined on same Domain");
774           data[i]=d.getDataC();           data[i]=d.getDataC();
# Line 844  void MeshAdapter::saveVTK(const std::str Line 782  void MeshAdapter::saveVTK(const std::str
782              strcpy(c_names[i],n.c_str());              strcpy(c_names[i],n.c_str());
783           }           }
784      }      }
785    #ifndef PASO_MPI    
786      Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);      Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
787      checkFinleyError();  #else
788        Finley_Mesh_saveVTK_MPIO(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
789    #endif
790    
791    checkFinleyError();
792      /* win32 refactor */
793      TMPMEMFREE(c_names);
794      TMPMEMFREE(data);
795      TMPMEMFREE(ptr_data);
796      for(int i=0;i<num_data;i++)
797      {
798        TMPMEMFREE(names[i]);
799      }
800      TMPMEMFREE(names);
801    
802      return;      return;
803  }  }
804                                                                                                                                                                                                                                                                                                                                                    
# Line 1041  int MeshAdapter::getSystemMatrixTypeId(c Line 994  int MeshAdapter::getSystemMatrixTypeId(c
994     return out;     return out;
995  }  }
996    
997  Data MeshAdapter::getX() const  escript::Data MeshAdapter::getX() const
998  {  {
999    return continuousFunction(asAbstractContinuousDomain()).getX();    return continuousFunction(asAbstractContinuousDomain()).getX();
1000  }  }
1001    
1002  Data MeshAdapter::getNormal() const  escript::Data MeshAdapter::getNormal() const
1003  {  {
1004    return functionOnBoundary(asAbstractContinuousDomain()).getNormal();    return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
1005  }  }
1006    
1007  Data MeshAdapter::getSize() const  escript::Data MeshAdapter::getSize() const
1008  {  {
1009    return function(asAbstractContinuousDomain()).getSize();    return function(asAbstractContinuousDomain()).getSize();
1010  }  }
1011    
1012  int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const  int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
1013  {  {
1014    int* tagList;    int out=0;
1015    int numTags;    Finley_Mesh* mesh=m_finleyMesh.get();
1016    getTagList(functionSpaceType, &tagList, &numTags);    switch (functionSpaceType) {
1017    return tagList[sampleNo];    case(Nodes):
1018        out=mesh->Nodes->Tag[sampleNo];
1019        break;
1020      case(Elements):
1021        out=mesh->Elements->Tag[sampleNo];
1022        break;
1023      case(FaceElements):
1024        out=mesh->FaceElements->Tag[sampleNo];
1025        break;
1026      case(Points):
1027        out=mesh->Points->Tag[sampleNo];
1028        break;
1029      case(ContactElementsZero):
1030        out=mesh->ContactElements->Tag[sampleNo];
1031        break;
1032      case(ContactElementsOne):
1033        out=mesh->ContactElements->Tag[sampleNo];
1034        break;
1035      case(DegreesOfFreedom):
1036        break;
1037      case(ReducedDegreesOfFreedom):
1038        break;
1039      default:
1040        stringstream temp;
1041        temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1042        throw FinleyAdapterException(temp.str());
1043        break;
1044      }
1045      return out;
1046    }
1047    int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1048    {
1049      int *out=0,i;
1050      Finley_Mesh* mesh=m_finleyMesh.get();
1051      switch (functionSpaceType) {
1052      case(Nodes):
1053        if (mesh->Nodes!=NULL) {
1054          out=mesh->Nodes->Id;
1055          break;
1056        }
1057      case(Elements):
1058        out=mesh->Elements->Id;
1059        break;
1060      case(FaceElements):
1061        out=mesh->FaceElements->Id;
1062        break;
1063      case(Points):
1064        out=mesh->Points->Id;
1065        break;
1066      case(ContactElementsZero):
1067        out=mesh->ContactElements->Id;
1068        break;
1069      case(ContactElementsOne):
1070        out=mesh->ContactElements->Id;
1071        break;
1072      case(DegreesOfFreedom):
1073        out=mesh->Nodes->degreeOfFreedomId;
1074        break;
1075      case(ReducedDegreesOfFreedom):
1076        out=mesh->Nodes->reducedDegreeOfFreedomId;
1077        break;
1078      default:
1079        stringstream temp;
1080        temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
1081        throw FinleyAdapterException(temp.str());
1082        break;
1083      }
1084      return out;
1085  }  }
1086    
1087  int MeshAdapter::getReferenceNoFromSampleNo(int functionSpaceType, int sampleNo) const  void MeshAdapter::setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const
1088  {  {
1089    int* referenceNoList;    Finley_Mesh* mesh=m_finleyMesh.get();
1090    int numReferenceNo;    escriptDataC tmp=mask.getDataC();
1091    getReferenceNoList(functionSpaceType, &referenceNoList, &numReferenceNo);    switch(functionSpaceType) {
1092    return referenceNoList[sampleNo];         case(Nodes):
1093              Finley_NodeFile_setTags(mesh->Nodes,newTag,&tmp);
1094              break;
1095           case(DegreesOfFreedom):
1096              throw FinleyAdapterException("Error - DegreesOfFreedom does not support tags");
1097              break;
1098           case(ReducedDegreesOfFreedom):
1099              throw FinleyAdapterException("Error - ReducedDegreesOfFreedom does not support tags");
1100              break;
1101           case(Elements):
1102              Finley_ElementFile_setTags(mesh->Elements,newTag,&tmp);
1103              break;
1104           case(FaceElements):
1105              Finley_ElementFile_setTags(mesh->FaceElements,newTag,&tmp);
1106              break;
1107           case(Points):
1108              Finley_ElementFile_setTags(mesh->Points,newTag,&tmp);
1109              break;
1110           case(ContactElementsZero):
1111              Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1112              break;
1113           case(ContactElementsOne):
1114              Finley_ElementFile_setTags(mesh->ContactElements,newTag,&tmp);
1115              break;
1116           default:
1117              stringstream temp;
1118              temp << "Error - Finley does not know anything about function space type " << functionSpaceType;
1119              throw FinleyAdapterException(temp.str());
1120      }
1121      checkFinleyError();
1122      return;
1123  }  }
1124    
1125    
1126  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.472  
changed lines
  Added in v.967

  ViewVC Help
Powered by ViewVC 1.1.26