/[escript]/branches/domexper/dudley/src/CPPAdapter/MeshAdapter.cpp
ViewVC logotype

Diff of /branches/domexper/dudley/src/CPPAdapter/MeshAdapter.cpp

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

revision 2385 by gross, Wed Apr 15 03:52:38 2009 UTC revision 2640 by jfenwick, Mon Aug 31 06:22:10 2009 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*******************************************************
3  *  *
4  * Copyright (c) 2003-2008 by University of Queensland  * Copyright (c) 2003-2009 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 96  bool MeshAdapter::onMasterProcessor() co Line 96  bool MeshAdapter::onMasterProcessor() co
96  }  }
97    
98    
99    #ifdef PASO_MPI
100      MPI_Comm
101    #else
102      unsigned int
103    #endif
104    MeshAdapter::getMPIComm() const
105    {
106    #ifdef PASO_MPI
107        return mesh->MPIInfo->comm;
108    #else
109        return 0;
110    #endif
111    }
112    
113    
114  Finley_Mesh* MeshAdapter::getFinley_Mesh() const {  Finley_Mesh* MeshAdapter::getFinley_Mesh() const {
115     return m_finleyMesh.get();     return m_finleyMesh.get();
116  }  }
# Line 109  void MeshAdapter::write(const std::strin Line 124  void MeshAdapter::write(const std::strin
124     TMPMEMFREE(fName);     TMPMEMFREE(fName);
125  }  }
126    
127  void MeshAdapter::Print_Mesh_Info(const bool full=false) const  void MeshAdapter::Print_Mesh_Info(const bool full) const
128  {  {
129     Finley_PrintMesh_Info(m_finleyMesh.get(), full);     Finley_PrintMesh_Info(m_finleyMesh.get(), full);
130  }  }
# Line 1453  void MeshAdapter::setNewX(const escript: Line 1468  void MeshAdapter::setNewX(const escript:
1468     const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(*(new_x.getFunctionSpace().getDomain()));     const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(*(new_x.getFunctionSpace().getDomain()));
1469     if (newDomain!=*this)     if (newDomain!=*this)
1470        throw FinleyAdapterException("Error - Illegal domain of new point locations");        throw FinleyAdapterException("Error - Illegal domain of new point locations");
1471     tmp = new_x.getDataC();     if ( new_x.getFunctionSpace() == continuousFunction(asAbstractContinuousDomain()) ) {
1472     Finley_Mesh_setCoordinates(mesh,&tmp);         tmp = new_x.getDataC();
1473           Finley_Mesh_setCoordinates(mesh,&tmp);
1474       } else {
1475           escript::Data new_x_inter=escript::Data( new_x,  continuousFunction(asAbstractContinuousDomain()) );
1476           tmp = new_x_inter.getDataC();
1477           Finley_Mesh_setCoordinates(mesh,&tmp);
1478       }
1479     checkFinleyError();     checkFinleyError();
1480  }  }
1481    
# Line 1512  void MeshAdapter::saveDX(const std::stri Line 1533  void MeshAdapter::saveDX(const std::stri
1533  //  //
1534  // saves mesh and optionally data arrays in VTK format  // saves mesh and optionally data arrays in VTK format
1535  //  //
1536  void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const  void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg,  const std::string& metadata, const std::string& metadata_schema) const
1537  {  {
1538     int num_data;     int num_data;
1539     char **names;     char **names;
# Line 1520  void MeshAdapter::saveVTK(const std::str Line 1541  void MeshAdapter::saveVTK(const std::str
1541     escriptDataC **ptr_data;     escriptDataC **ptr_data;
1542    
1543     extractArgsFromDict(arg, num_data, names, data, ptr_data);     extractArgsFromDict(arg, num_data, names, data, ptr_data);
1544     Finley_Mesh_saveVTK(filename.c_str(), m_finleyMesh.get(), num_data, names, ptr_data);     Finley_Mesh_saveVTK(filename.c_str(), m_finleyMesh.get(), num_data, names, ptr_data, metadata.c_str(), metadata_schema.c_str());
1545     checkFinleyError();     checkFinleyError();
1546    
1547     /* win32 refactor */     /* win32 refactor */
# Line 1579  SystemMatrixAdapter MeshAdapter::newSyst Line 1600  SystemMatrixAdapter MeshAdapter::newSyst
1600  #endif  #endif
1601     }     }
1602     else {     else {
1603        fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);        fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize,FALSE);
1604     }     }
1605     checkPasoError();     checkPasoError();
1606     Paso_SystemMatrixPattern_free(fsystemMatrixPattern);     Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
# Line 1655  bool MeshAdapter::isCellOriented(int fun Line 1676  bool MeshAdapter::isCellOriented(int fun
1676     return false;     return false;
1677  }  }
1678    
1679    bool
1680    MeshAdapter::commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const
1681    {
1682       /* The idea is to use equivalence classes. [Types which can be interpolated back and forth]
1683        class 1: DOF <-> Nodes
1684        class 2: ReducedDOF <-> ReducedNodes
1685        class 3: Points
1686        class 4: Elements
1687        class 5: ReducedElements
1688        class 6: FaceElements
1689        class 7: ReducedFaceElements
1690        class 8: ContactElementZero <-> ContactElementOne
1691        class 9: ReducedContactElementZero <-> ReducedContactElementOne
1692    
1693       There is also a set of lines. Interpolation is possible down a line but not between lines.
1694       class 1 and 2 belong to all lines so aren't considered.
1695        line 0: class 3
1696        line 1: class 4,5
1697        line 2: class 6,7
1698        line 3: class 8,9
1699    
1700       For classes with multiple members (eg class 2) we have vars to record if there is at least one instance.
1701       eg hasnodes is true if we have at least one instance of Nodes.
1702       */
1703        if (fs.size()==0)
1704        {
1705        return false;
1706        }
1707        std::vector<int> hasclass(10);
1708        std::vector<int> hasline(4);    
1709        bool hasnodes=false;
1710        bool hasrednodes=false;
1711        bool hascez=false;
1712        bool hasrcez=false;
1713        for (int i=0;i<fs.size();++i)
1714        {
1715        switch(fs[i])
1716        {
1717        case(Nodes):   hasnodes=true;   // no break is deliberate
1718        case(DegreesOfFreedom):
1719            hasclass[1]=1;
1720            break;
1721        case(ReducedNodes):    hasrednodes=true;    // no break is deliberate
1722        case(ReducedDegreesOfFreedom):
1723            hasclass[2]=1;
1724            break;
1725        case(Points):
1726            hasline[0]=1;
1727            hasclass[3]=1;
1728            break;
1729        case(Elements):
1730            hasclass[4]=1;
1731            hasline[1]=1;
1732            break;
1733        case(ReducedElements):
1734            hasclass[5]=1;
1735            hasline[1]=1;
1736            break;
1737        case(FaceElements):
1738            hasclass[6]=1;
1739            hasline[2]=1;
1740            break;
1741        case(ReducedFaceElements):
1742            hasclass[7]=1;
1743            hasline[2]=1;
1744            break;
1745        case(ContactElementsZero):  hascez=true;    // no break is deliberate
1746        case(ContactElementsOne):
1747            hasclass[8]=1;
1748            hasline[3]=1;
1749            break;
1750        case(ReducedContactElementsZero):   hasrcez=true;   // no break is deliberate
1751        case(ReducedContactElementsOne):
1752            hasclass[9]=1;
1753            hasline[3]=1;
1754            break;
1755        default:
1756            return false;
1757        }
1758        }
1759        int totlines=hasline[0]+hasline[1]+hasline[2]+hasline[3];
1760        // fail if we have more than one leaf group
1761    
1762        if (totlines>1)
1763        {
1764        return false;   // there are at least two branches we can't interpolate between
1765        }
1766        else if (totlines==1)
1767        {
1768        if (hasline[0]==1)      // we have points
1769        {
1770            resultcode=Points;
1771        }
1772        else if (hasline[1]==1)
1773        {
1774            if (hasclass[5]==1)
1775            {
1776            resultcode=ReducedElements;
1777            }
1778            else
1779            {
1780            resultcode=Elements;
1781            }
1782        }
1783        else if (hasline[2]==1)
1784        {
1785            if (hasclass[7]==ReducedFaceElements)
1786            {
1787            resultcode=ReducedFaceElements;
1788            }
1789            else
1790            {
1791            resultcode=FaceElements;
1792            }
1793        }
1794        else    // so we must be in line3
1795        {
1796            if (hasclass[9]==1)
1797            {
1798            // need something from class 9
1799            resultcode=(hasrcez?ReducedContactElementsZero:ReducedContactElementsOne);
1800            }
1801            else
1802            {
1803            // something from class 8
1804            resultcode=(hascez?ContactElementsZero:ContactElementsOne);
1805            }
1806        }
1807        }
1808        else    // totlines==0
1809        {
1810        if (hasclass[2]==1)
1811        {
1812            // something from class 2
1813            resultcode=(hasrednodes?ReducedNodes:ReducedDegreesOfFreedom);
1814        }
1815        else
1816        {   // something from class 1
1817            resultcode=(hasnodes?Nodes:DegreesOfFreedom);
1818        }
1819        }
1820        return true;
1821    }
1822    
1823  bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const  bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const
1824  {  {
1825     switch(functionSpaceType_source) {     switch(functionSpaceType_source) {
1826     case(Nodes):     case(Nodes):
1827     switch(functionSpaceType_target) {      switch(functionSpaceType_target) {
1828     case(Nodes):      case(Nodes):
1829     case(ReducedNodes):      case(ReducedNodes):
1830     case(ReducedDegreesOfFreedom):      case(ReducedDegreesOfFreedom):
1831     case(DegreesOfFreedom):      case(DegreesOfFreedom):
1832     case(Elements):      case(Elements):
1833     case(ReducedElements):      case(ReducedElements):
1834     case(FaceElements):      case(FaceElements):
1835     case(ReducedFaceElements):      case(ReducedFaceElements):
1836     case(Points):      case(Points):
1837     case(ContactElementsZero):      case(ContactElementsZero):
1838     case(ReducedContactElementsZero):      case(ReducedContactElementsZero):
1839     case(ContactElementsOne):      case(ContactElementsOne):
1840     case(ReducedContactElementsOne):      case(ReducedContactElementsOne):
1841     return true;      return true;
1842     default:      default:
1843        stringstream temp;            stringstream temp;
1844        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;            temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1845        throw FinleyAdapterException(temp.str());            throw FinleyAdapterException(temp.str());
1846     }     }
1847     break;     break;
1848     case(ReducedNodes):     case(ReducedNodes):
1849     switch(functionSpaceType_target) {      switch(functionSpaceType_target) {
1850     case(ReducedNodes):      case(ReducedNodes):
1851     case(ReducedDegreesOfFreedom):      case(ReducedDegreesOfFreedom):
1852     case(Elements):      case(Elements):
1853     case(ReducedElements):      case(ReducedElements):
1854     case(FaceElements):      case(FaceElements):
1855     case(ReducedFaceElements):      case(ReducedFaceElements):
1856     case(Points):      case(Points):
1857     case(ContactElementsZero):      case(ContactElementsZero):
1858     case(ReducedContactElementsZero):      case(ReducedContactElementsZero):
1859     case(ContactElementsOne):      case(ContactElementsOne):
1860     case(ReducedContactElementsOne):      case(ReducedContactElementsOne):
1861     return true;      return true;
1862     case(Nodes):      case(Nodes):
1863     case(DegreesOfFreedom):      case(DegreesOfFreedom):
1864     return false;      return false;
1865     default:      default:
1866        stringstream temp;          stringstream temp;
1867        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;          temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1868        throw FinleyAdapterException(temp.str());          throw FinleyAdapterException(temp.str());
1869     }     }
1870     break;     break;
1871     case(Elements):     case(Elements):
1872     if (functionSpaceType_target==Elements) {      if (functionSpaceType_target==Elements) {
1873        return true;        return true;
1874     } else if (functionSpaceType_target==ReducedElements) {      } else if (functionSpaceType_target==ReducedElements) {
1875        return true;        return true;
1876     } else {          } else {
1877        return false;            return false;
1878     }          }
1879     case(ReducedElements):     case(ReducedElements):
1880     if (functionSpaceType_target==ReducedElements) {      if (functionSpaceType_target==ReducedElements) {
1881        return true;        return true;
1882     } else {      } else {
1883        return false;            return false;
1884     }      }
1885     case(FaceElements):     case(FaceElements):
1886     if (functionSpaceType_target==FaceElements) {      if (functionSpaceType_target==FaceElements) {
1887        return true;              return true;
1888     } else if (functionSpaceType_target==ReducedFaceElements) {      } else if (functionSpaceType_target==ReducedFaceElements) {
1889        return true;              return true;
1890     } else {      } else {
1891        return false;              return false;
1892     }      }
1893     case(ReducedFaceElements):     case(ReducedFaceElements):
1894     if (functionSpaceType_target==ReducedFaceElements) {      if (functionSpaceType_target==ReducedFaceElements) {
1895        return true;              return true;
1896     } else {      } else {
1897        return false;          return false;
1898     }      }
1899     case(Points):     case(Points):
1900     if (functionSpaceType_target==Points) {      if (functionSpaceType_target==Points) {
1901        return true;              return true;
1902     } else {      } else {
1903        return false;              return false;
1904     }      }
1905     case(ContactElementsZero):     case(ContactElementsZero):
1906     case(ContactElementsOne):     case(ContactElementsOne):
1907     if (functionSpaceType_target==ContactElementsZero || functionSpaceType_target==ContactElementsOne) {      if (functionSpaceType_target==ContactElementsZero || functionSpaceType_target==ContactElementsOne) {
1908        return true;              return true;
1909     } else if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {      } else if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1910        return true;              return true;
1911     } else {      } else {
1912        return false;              return false;
1913     }      }
1914     case(ReducedContactElementsZero):     case(ReducedContactElementsZero):
1915     case(ReducedContactElementsOne):     case(ReducedContactElementsOne):
1916     if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {      if (functionSpaceType_target==ReducedContactElementsZero || functionSpaceType_target==ReducedContactElementsOne) {
1917        return true;              return true;
1918     } else {      } else {
1919        return false;              return false;
1920     }      }
1921     case(DegreesOfFreedom):     case(DegreesOfFreedom):
1922     switch(functionSpaceType_target) {      switch(functionSpaceType_target) {
1923     case(ReducedDegreesOfFreedom):      case(ReducedDegreesOfFreedom):
1924     case(DegreesOfFreedom):      case(DegreesOfFreedom):
1925     case(Nodes):      case(Nodes):
1926     case(ReducedNodes):      case(ReducedNodes):
1927     case(Elements):      case(Elements):
1928     case(ReducedElements):      case(ReducedElements):
1929     case(Points):      case(Points):
1930     case(FaceElements):      case(FaceElements):
1931     case(ReducedFaceElements):      case(ReducedFaceElements):
1932     case(ContactElementsZero):      case(ContactElementsZero):
1933     case(ReducedContactElementsZero):      case(ReducedContactElementsZero):
1934     case(ContactElementsOne):      case(ContactElementsOne):
1935     case(ReducedContactElementsOne):      case(ReducedContactElementsOne):
1936     return true;      return true;
1937     default:      default:
1938        stringstream temp;          stringstream temp;
1939        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;          temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1940        throw FinleyAdapterException(temp.str());          throw FinleyAdapterException(temp.str());
1941     }      }
1942     break;      break;
1943     case(ReducedDegreesOfFreedom):     case(ReducedDegreesOfFreedom):
1944     switch(functionSpaceType_target) {     switch(functionSpaceType_target) {
1945     case(ReducedDegreesOfFreedom):      case(ReducedDegreesOfFreedom):
1946     case(ReducedNodes):      case(ReducedNodes):
1947     case(Elements):      case(Elements):
1948     case(ReducedElements):      case(ReducedElements):
1949     case(FaceElements):      case(FaceElements):
1950     case(ReducedFaceElements):      case(ReducedFaceElements):
1951     case(Points):      case(Points):
1952     case(ContactElementsZero):      case(ContactElementsZero):
1953     case(ReducedContactElementsZero):      case(ReducedContactElementsZero):
1954     case(ContactElementsOne):      case(ContactElementsOne):
1955     case(ReducedContactElementsOne):      case(ReducedContactElementsOne):
1956     return true;      return true;
1957     case(Nodes):      case(Nodes):
1958     case(DegreesOfFreedom):      case(DegreesOfFreedom):
1959     return false;      return false;
1960     default:      default:
1961        stringstream temp;          stringstream temp;
1962        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;          temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
1963        throw FinleyAdapterException(temp.str());          throw FinleyAdapterException(temp.str());
1964     }      }
1965     break;      break;
1966     default:     default:
1967        stringstream temp;        stringstream temp;
1968        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;        temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;
# Line 1858  escript::Data MeshAdapter::getSize() con Line 2023  escript::Data MeshAdapter::getSize() con
2023     return escript::function(asAbstractContinuousDomain()).getSize();     return escript::function(asAbstractContinuousDomain()).getSize();
2024  }  }
2025    
2026  int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const  const int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
2027  {  {
2028     int *out = NULL;     int *out = NULL;
2029     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
# Line 2095  int MeshAdapter::getNumberOfTagsInUse(in Line 2260  int MeshAdapter::getNumberOfTagsInUse(in
2260    }    }
2261    return numTags;    return numTags;
2262  }  }
2263  int* MeshAdapter::borrowListOfTagsInUse(int functionSpaceCode) const  
2264    const int* MeshAdapter::borrowListOfTagsInUse(int functionSpaceCode) const
2265  {  {
2266    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
2267    index_t* tags=NULL;    index_t* tags=NULL;
# Line 2161  bool MeshAdapter::canTag(int functionSpa Line 2327  bool MeshAdapter::canTag(int functionSpa
2327    }    }
2328  }  }
2329    
2330    AbstractDomain::StatusType MeshAdapter::getStatus() const
2331    {
2332      Finley_Mesh* mesh=m_finleyMesh.get();
2333      return Finley_Mesh_getStatus(mesh);
2334    }
2335    
2336    
2337    
2338  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.2385  
changed lines
  Added in v.2640

  ViewVC Help
Powered by ViewVC 1.1.26