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

Diff of /branches/diaplayground/finley/src/CPPAdapter/MeshAdapter.cpp

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

revision 3546 by jfenwick, Thu Jul 28 01:32:20 2011 UTC revision 3911 by jfenwick, Thu Jun 14 01:01:03 2012 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*******************************************************
3  *  *
4  * Copyright (c) 2003-2010 by University of Queensland  * Copyright (c) 2003-2012 by University of Queensland
5  * Earth Systems Science Computational Center (ESSCC)  * Earth Systems Science Computational Center (ESSCC)
6  * http://www.uq.edu.au/esscc  * http://www.uq.edu.au/esscc
7  *  *
# Line 11  Line 11 
11  *  *
12  *******************************************************/  *******************************************************/
13    
14    #include <pasowrap/PasoException.h>
15    #include <pasowrap/TransportProblemAdapter.h>
16  #include "MeshAdapter.h"  #include "MeshAdapter.h"
17  #include "escript/Data.h"  #include "escript/Data.h"
18  #include "escript/DataFactory.h"  #include "escript/DataFactory.h"
# Line 26  extern "C" { Line 27  extern "C" {
27    
28  using namespace std;  using namespace std;
29  using namespace escript;  using namespace escript;
30    using namespace paso;
31    
32  namespace finley {  namespace finley {
33    
# Line 777  void MeshAdapter::addPDEToSystem( Line 779  void MeshAdapter::addPDEToSystem(
779     SystemMatrixAdapter* smat=dynamic_cast<SystemMatrixAdapter*>(&mat);     SystemMatrixAdapter* smat=dynamic_cast<SystemMatrixAdapter*>(&mat);
780     if (smat==0)     if (smat==0)
781     {     {
782      throw FinleyAdapterException("finley only supports its own system matrices.");      throw FinleyAdapterException("finley only supports Paso system matrices.");
783     }     }
784     escriptDataC _rhs=rhs.getDataC();     escriptDataC _rhs=rhs.getDataC();
785     escriptDataC _A  =A.getDataC();     escriptDataC _A  =A.getDataC();
# Line 875  void MeshAdapter::addPDEToTransportProbl Line 877  void MeshAdapter::addPDEToTransportProbl
877     TransportProblemAdapter* tpa=dynamic_cast<TransportProblemAdapter*>(&tp);     TransportProblemAdapter* tpa=dynamic_cast<TransportProblemAdapter*>(&tp);
878     if (tpa==0)     if (tpa==0)
879     {     {
880      throw FinleyAdapterException("finley only supports its own transport problems.");      throw FinleyAdapterException("finley only supports Paso transport problems.");
881     }     }
882    
883    
# Line 1663  ASM_ptr MeshAdapter::newSystemMatrix( Line 1665  ASM_ptr MeshAdapter::newSystemMatrix(
1665  // creates a TransportProblemAdapter  // creates a TransportProblemAdapter
1666  //  //
1667  ATP_ptr MeshAdapter::newTransportProblem(  ATP_ptr MeshAdapter::newTransportProblem(
                                                          const bool useBackwardEuler,  
1668                                                           const int blocksize,                                                           const int blocksize,
1669                                                           const escript::FunctionSpace& functionspace,                                                           const escript::FunctionSpace& functionspace,
1670                                                           const int type) const                                                           const int type) const
# Line 1686  ATP_ptr MeshAdapter::newTransportProblem Line 1687  ATP_ptr MeshAdapter::newTransportProblem
1687     Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);     Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);
1688     checkFinleyError();     checkFinleyError();
1689     Paso_TransportProblem* transportProblem;     Paso_TransportProblem* transportProblem;
1690     transportProblem=Paso_TransportProblem_alloc(useBackwardEuler,fsystemMatrixPattern,blocksize);     transportProblem=Paso_TransportProblem_alloc(fsystemMatrixPattern,blocksize);
1691     checkPasoError();     checkPasoError();
1692     Paso_SystemMatrixPattern_free(fsystemMatrixPattern);     Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1693     TransportProblemAdapter* tpa=new TransportProblemAdapter(transportProblem,useBackwardEuler,blocksize,functionspace);     TransportProblemAdapter* tpa=new TransportProblemAdapter(transportProblem,blocksize,functionspace);
1694     return ATP_ptr(tpa);     return ATP_ptr(tpa);
1695  //   return TransportProblemAdapter(transportProblem,useBackwardEuler,blocksize,functionspace);  //   return TransportProblemAdapter(transportProblem, blocksize,functionspace);
1696  }  }
1697    
1698  //  //
# Line 2050  bool MeshAdapter::operator!=(const Abstr Line 2051  bool MeshAdapter::operator!=(const Abstr
2051  int MeshAdapter::getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const  int MeshAdapter::getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const
2052  {  {
2053     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
2054     int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(preconditioner), SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0, mesh->MPIInfo);     return SystemMatrixAdapter::getSystemMatrixTypeId(solver, preconditioner,
2055     checkPasoError();             package, symmetry, mesh->MPIInfo);
    return out;  
2056  }  }
2057    
2058  int MeshAdapter::getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const  int MeshAdapter::getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const
2059  {  {
2060     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
2061     int out=Paso_TransportProblem_getTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(preconditioner), SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0, mesh->MPIInfo);     return TransportProblemAdapter::getTransportTypeId(solver, preconditioner,
2062     checkPasoError();             package, symmetry, mesh->MPIInfo);
    return out;  
2063  }  }
2064    
2065  escript::Data MeshAdapter::getX() const  escript::Data MeshAdapter::getX() const
# Line 2240  void MeshAdapter::setTagMap(const string Line 2240  void MeshAdapter::setTagMap(const string
2240  {  {
2241     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
2242     Finley_Mesh_addTagMap(mesh, name.c_str(),tag);     Finley_Mesh_addTagMap(mesh, name.c_str(),tag);
2243     checkPasoError();     checkFinleyError();
2244     // throwStandardException("MeshAdapter::set TagMap is not implemented.");     // throwStandardException("MeshAdapter::set TagMap is not implemented.");
2245  }  }
2246    
# Line 2249  int MeshAdapter::getTag(const string& na Line 2249  int MeshAdapter::getTag(const string& na
2249     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
2250     int tag=0;     int tag=0;
2251     tag=Finley_Mesh_getTag(mesh, name.c_str());     tag=Finley_Mesh_getTag(mesh, name.c_str());
2252     checkPasoError();     checkFinleyError();
2253     // throwStandardException("MeshAdapter::getTag is not implemented.");     // throwStandardException("MeshAdapter::getTag is not implemented.");
2254     return tag;     return tag;
2255  }  }
# Line 2437  ReferenceElementSetWrapper::~ReferenceEl Line 2437  ReferenceElementSetWrapper::~ReferenceEl
2437    Finley_ReferenceElementSet_dealloc(m_refSet);    Finley_ReferenceElementSet_dealloc(m_refSet);
2438  }  }
2439    
2440  void MeshAdapter:: addDiracPoints(const boost::python::list& points, const boost::python::list& tags) const  // points will be flattened
2441    void MeshAdapter:: addDiracPoints(const std::vector<double>& points, const std::vector<int>& tags) const
2442  {  {
2443        const int dim = getDim();        const int dim = getDim();
2444        int numPoints=boost::python::extract<int>(points.attr("__len__")());        int numPoints=points.size()/dim;
2445        int numTags=boost::python::extract<int>(tags.attr("__len__")());        int numTags=tags.size();
2446        Finley_Mesh* mesh=m_finleyMesh.get();        Finley_Mesh* mesh=m_finleyMesh.get();
2447                
2448          if ( points.size() % dim != 0 )
2449          {
2450        throw FinleyAdapterException("Error - number of coords does not appear to be a multiple of dimension.");
2451          }
2452          
2453        if  ( (numTags > 0) && ( numPoints !=  numTags ) )        if  ( (numTags > 0) && ( numPoints !=  numTags ) )
2454       throw FinleyAdapterException("Error - if tags are given number of tags and points must match.");       throw FinleyAdapterException("Error - if tags are given number of tags and points must match.");
2455                
# Line 2451  void MeshAdapter:: addDiracPoints(const Line 2457  void MeshAdapter:: addDiracPoints(const
2457        int*    tags_ptr= TMPMEMALLOC(numPoints, int);        int*    tags_ptr= TMPMEMALLOC(numPoints, int);
2458                
2459        for (int i=0;i<numPoints;++i) {        for (int i=0;i<numPoints;++i) {
2460         int tag_id=-1;         points_ptr[ i * dim     ] = points[i * dim ];
2461         int numComps=boost::python::extract<int>(points[i].attr("__len__")());         if ( dim > 1 ) points_ptr[ i * dim + 1 ] = points[i * dim + 1];
2462         if  ( numComps !=   dim ) {         if ( dim > 2 ) points_ptr[ i * dim + 2 ] = points[i * dim + 2];    
2463                 stringstream temp;                       tags_ptr[i]=tags[i];
                temp << "Error - illegal number of components " << numComps << " for point " << i;                
                throw FinleyAdapterException(temp.str());  
        }  
        points_ptr[ i * dim     ] = boost::python::extract<double>(points[i][0]);  
        if ( dim > 1 ) points_ptr[ i * dim + 1 ] = boost::python::extract<double>(points[i][1]);  
        if ( dim > 2 ) points_ptr[ i * dim + 2 ] = boost::python::extract<double>(points[i][2]);  
         
        if ( numTags > 0) {  
               boost::python::extract<string> ex_str(tags[i]);  
           if  ( ex_str.check() ) {  
               tag_id=getTag( ex_str());  
           } else {  
                boost::python::extract<int> ex_int(tags[i]);  
                if ( ex_int.check() ) {  
                    tag_id=ex_int();  
                } else {  
                 stringstream temp;            
                     temp << "Error - unable to extract tag for point " << i;  
                 throw FinleyAdapterException(temp.str());  
               }  
           }  
        }        
            tags_ptr[i]=tag_id;  
2464        }        }
2465                
2466        Finley_Mesh_addPoints(mesh, numPoints, points_ptr, tags_ptr);        Finley_Mesh_addPoints(mesh, numPoints, points_ptr, tags_ptr);
2467        checkPasoError();        checkFinleyError();
2468                
2469        TMPMEMFREE(points_ptr);        TMPMEMFREE(points_ptr);
2470        TMPMEMFREE(tags_ptr);        TMPMEMFREE(tags_ptr);
2471  }  }
2472    
2473    
2474    // void MeshAdapter:: addDiracPoints(const boost::python::list& points, const boost::python::list& tags) const
2475    // {
2476    //       const int dim = getDim();
2477    //       int numPoints=boost::python::extract<int>(points.attr("__len__")());
2478    //       int numTags=boost::python::extract<int>(tags.attr("__len__")());
2479    //       Finley_Mesh* mesh=m_finleyMesh.get();
2480    //      
2481    //       if  ( (numTags > 0) && ( numPoints !=  numTags ) )
2482    //   throw FinleyAdapterException("Error - if tags are given number of tags and points must match.");
2483    //      
2484    //       double* points_ptr=TMPMEMALLOC(numPoints * dim, double);
2485    //       int*    tags_ptr= TMPMEMALLOC(numPoints, int);
2486    //      
2487    //       for (int i=0;i<numPoints;++i) {
2488    //     int tag_id=-1;
2489    //     int numComps=boost::python::extract<int>(points[i].attr("__len__")());
2490    //     if  ( numComps !=   dim ) {
2491    //                stringstream temp;            
2492    //                temp << "Error - illegal number of components " << numComps << " for point " << i;              
2493    //                throw FinleyAdapterException(temp.str());
2494    //     }
2495    //     points_ptr[ i * dim     ] = boost::python::extract<double>(points[i][0]);
2496    //     if ( dim > 1 ) points_ptr[ i * dim + 1 ] = boost::python::extract<double>(points[i][1]);
2497    //     if ( dim > 2 ) points_ptr[ i * dim + 2 ] = boost::python::extract<double>(points[i][2]);
2498    //    
2499    //     if ( numTags > 0) {
2500    //            boost::python::extract<string> ex_str(tags[i]);
2501    //        if  ( ex_str.check() ) {
2502    //            tag_id=getTag( ex_str());
2503    //        } else {
2504    //             boost::python::extract<int> ex_int(tags[i]);
2505    //             if ( ex_int.check() ) {
2506    //                 tag_id=ex_int();
2507    //             } else {
2508    //              stringstream temp;          
2509    //                  temp << "Error - unable to extract tag for point " << i;
2510    //              throw FinleyAdapterException(temp.str());
2511    //            }
2512    //        }
2513    //     }      
2514    //            tags_ptr[i]=tag_id;
2515    //       }
2516    //      
2517    //       Finley_Mesh_addPoints(mesh, numPoints, points_ptr, tags_ptr);
2518    //       checkPasoError();
2519    //      
2520    //       TMPMEMFREE(points_ptr);
2521    //       TMPMEMFREE(tags_ptr);
2522    // }
2523    
2524    /*
2525  void MeshAdapter:: addDiracPoint( const boost::python::list& point, const int tag) const  void MeshAdapter:: addDiracPoint( const boost::python::list& point, const int tag) const
2526  {    {  
2527      boost::python::list points =  boost::python::list();      boost::python::list points =  boost::python::list();
# Line 2495  void MeshAdapter:: addDiracPoint( const Line 2530  void MeshAdapter:: addDiracPoint( const
2530      tags.append(tag);      tags.append(tag);
2531      addDiracPoints(points, tags);      addDiracPoints(points, tags);
2532  }  }
2533    */
2534    
2535    /*
2536  void MeshAdapter:: addDiracPointWithTagName( const boost::python::list& point, const std::string& tag) const  void MeshAdapter:: addDiracPointWithTagName( const boost::python::list& point, const std::string& tag) const
2537  {  {
2538          boost::python::list points =   boost::python::list();          boost::python::list points =   boost::python::list();
# Line 2502  void MeshAdapter:: addDiracPointWithTagN Line 2540  void MeshAdapter:: addDiracPointWithTagN
2540          points.append(point);          points.append(point);
2541          tags.append(tag);          tags.append(tag);
2542          addDiracPoints(points, tags);          addDiracPoints(points, tags);
2543  }    }
2544    */
2545  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.3546  
changed lines
  Added in v.3911

  ViewVC Help
Powered by ViewVC 1.1.26