From: vbd Date: Thu, 6 Jan 2011 17:27:05 +0000 (+0000) Subject: addition of the MEDPartitioner tool (rewriting of medsplitter with MEDCoupling library) X-Git-Tag: MEDPartitioner_first_compilable_version~15 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=7ef298467d2619d8fd6a774c9f1bf8aa1c90ee5e;p=tools%2Fmedcoupling.git addition of the MEDPartitioner tool (rewriting of medsplitter with MEDCoupling library) --- diff --git a/src/MEDPartitioner/MEDPARTITIONER.hxx b/src/MEDPartitioner/MEDPARTITIONER.hxx new file mode 100755 index 000000000..80e64b033 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER.hxx @@ -0,0 +1,38 @@ + // Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +// File : MEDPARTITIONER.hxx +// Author : Alexander A. BORODIN +// Module : MED +// Exporting/Importing defines for Windows Platform +// +#ifndef MEDPARTITIONER_HXX_ +#define MEDPARTITIONER_HXX_ + +#ifdef WNT +# if defined MEDPARTITIONER_EXPORTS || defined medsplitter_EXPORTS +# define MEDPARTITIONER_EXPORT __declspec( dllexport ) +# else +# define MEDPARTITIONER_EXPORT __declspec( dllimport ) +# endif +#else +# define MEDPARTITIONER_EXPORT +#endif + +#endif //MEDPARTITIONER_HXX_ diff --git a/src/MEDPartitioner/MEDPARTITIONER_Graph.cxx b/src/MEDPartitioner/MEDPARTITIONER_Graph.cxx new file mode 100644 index 000000000..607aef3be --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_Graph.cxx @@ -0,0 +1,47 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +#include "MEDPARTITIONER_Graph.hxx" + +using namespace MEDPARTITIONER; + +Graph::Graph(MEDPARTITIONER::MEDSKYLINEARRAY* array, int* edgeweight):m_graph(array),m_partition(0),m_edgeweight(edgeweight),m_cellweight(0) +{ +} + +//================================================================================ +/*! + * \brief Destructor + */ +//================================================================================ + +Graph::~Graph() +{ + if (m_partition) + { + delete m_partition; + m_partition=0; + } + if (m_graph) + { + delete m_graph; + m_graph=0; + } +} + diff --git a/src/MEDPartitioner/MEDPARTITIONER_Graph.hxx b/src/MEDPartitioner/MEDPARTITIONER_Graph.hxx new file mode 100644 index 000000000..50af471fa --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_Graph.hxx @@ -0,0 +1,67 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +#ifndef MEDPARTITIONER_GRAPH_HXX_ +#define MEDPARTITIONER_GRAPH_HXX_ + +#include "MEDPARTITIONER.hxx" +#include "MEDPARTITIONER_SkyLineArray.hxx" + +namespace MEDPARTITIONER { + +class ParaDomainSelector; +class MEDPARTITIONER_EXPORT Graph +{ + public: + typedef enum {METIS,SCOTCH} splitter_type; + + Graph(){}; + + //creates a graph from a SKYLINEARRAY + Graph(MEDPARTITIONER::MEDSKYLINEARRAY* graph, int* edgeweight=0); + + virtual ~Graph(); + + void setEdgesWeights(int* edgeweight){m_edgeweight=edgeweight;} + void setVerticesWeights(int* cellweight){m_cellweight=cellweight;} + + //computes partitioning of the graph + virtual void partGraph(int ndomain, const std::string&, ParaDomainSelector* sel=0)=0; + + //! returns the partitioning + const int* getPart() const {return m_partition->getValue();} + + //! returns the number of graph vertices (which can correspond to the cells in the mesh!) + int nbVertices() const {return m_graph->getNumberOf();} + + const MEDPARTITIONER::MEDSKYLINEARRAY* getGraph() const {return m_graph;} + + protected: + + MEDPARTITIONER::MEDSKYLINEARRAY* m_graph; + + MEDPARTITIONER::MEDSKYLINEARRAY* m_partition; + + int* m_edgeweight; + + int* m_cellweight; + }; + +} +#endif /*GRAPH_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollection.cxx b/src/MEDPartitioner/MEDPARTITIONER_MESHCollection.cxx new file mode 100644 index 000000000..65addd7cd --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollection.cxx @@ -0,0 +1,2159 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +#include "MEDMEM_Exception.hxx" + +#include "MEDCouplingUMesh.hxx" +#include "MEDPARTITIONER_utils.hxx" + +#include "MEDPARTITIONER_Graph.hxx" + +#include "MEDPARTITIONER_Topology.hxx" +#include "MEDPARTITIONER_ParallelTopology.hxx" +#include "MEDPARTITIONER_SequentialTopology.hxx" +#include "MEDPARTITIONER_ParaDomainSelector.hxx" +//#include "MEDPARTITIONER_MeshSendReceive.hxx" +//#include "MEDPARTITIONER_JointExchangeData.hxx" + +#include "MEDPARTITIONER_MESHCollection.hxx" +#include "MEDPARTITIONER_MESHCollectionDriver.hxx" +#include "MEDPARTITIONER_MESHCollectionMedXMLDriver.hxx" +#include "MEDPARTITIONER_MESHCollectionMedAsciiDriver.hxx" + +#include "MEDPARTITIONER_UserGraph.hxx" + +#ifdef ENABLE_METIS +#include "MEDPARTITIONER_METISGraph.hxx" +#endif +#ifdef ENABLE_SCOTCH +#include "MEDPARTITIONER_SCOTCHGraph.hxx" +#endif + +#include +#include + +#ifndef WNT +# include +#else +# include +#endif +#include + +#include +#include + +using namespace MEDPARTITIONER; + +#ifndef WNT +using namespace __gnu_cxx; +#else +using namespace std; +#endif + +//template inclusion +//#include "MEDPARTITIONER_MESHCollection.H" + +MESHCollection::MESHCollection() + : _topology(0), + _owns_topology(false), + _driver(0), + _domain_selector( 0 ), + _i_non_empty_mesh(-1), + _driver_type(MEDPARTITIONER::MedXML), + _subdomain_boundary_creates(false), + _family_splitting(false), + _create_empty_groups(false) +{ +} + +/*!constructor creating a new mesh collection (mesh series + topology) + *from an old collection and a new topology + * + * On output, the constructor has built the meshes corresponding to the new mesh collection. + * The new topology has been updated so that face and node mappings are included. + * The families have been cast to their projections in the new topology. + * + * \param initial_collection collection from which the data (coordinates, connectivity) are taken + * \param topology topology containing the cell mappings + */ + +MESHCollection::MESHCollection(MESHCollection& initial_collection, Topology* topology, bool family_splitting, bool create_empty_groups) + : _name(initial_collection._name), + _topology(topology), + _owns_topology(false), + _driver(0), + _domain_selector( initial_collection._domain_selector ), + _i_non_empty_mesh(-1), + _driver_type(MEDPARTITIONER::MedXML), + _subdomain_boundary_creates(false), + _family_splitting(family_splitting), + _create_empty_groups(create_empty_groups) +{ + string mesh_name = initial_collection.getName(); + _mesh.resize(_topology->nbDomain()); + + //splitting the initial domains into smaller bits + + std::vector > splitMeshes; + splitMeshes.resize(topology->nbDomain()); + for (int inew=0; inewnbDomain();inew++) + splitMeshes[inew].resize(initial_collection.getTopology()->nbDomain()); + + for (int iold=0; ioldnbDomain();iold++) + { + if (!isParallelMode() || initial_collection._domain_selector->isMyDomain(iold)) + { + int size=(initial_collection._mesh)[iold]->getNumberOfCells(); + std::vector globalids(size); + initial_collection.getTopology()->getCellList(iold, &globalids[0]); + std::vector ilocalnew(size); + std::vector ipnew(size); + topology->convertGlobalCellList(&globalids[0],size,&ilocalnew[0],&ipnew[0]); + std::vector > ids(topology->nbDomain()); + for (int i=0; inbDomain();inew++) + { + splitMeshes[inew][iold]=(ParaMEDMEM::MEDCouplingUMesh*)(initial_collection.getMesh())[iold]->buildPartOfMySelf(&ids[inew][0],&ids[inew][0]+ids[inew].size(),true); + // cout <<"small domain "<getNumberOfCells()<<" cells"<nbDomain() ;inew++) + { + std::vector meshes(splitMeshes[inew].size()); + for (int i=0; i< splitMeshes[inew].size();i++) + meshes[i]=splitMeshes[inew][i]; + + if (!isParallelMode()||_domain_selector->isMyDomain(inew)) + { + _mesh[inew]=ParaMEDMEM::MEDCouplingUMesh::mergeUMeshes(meshes); + _mesh[inew]->zipCoords(); + + // std::cout<<"creating mesh #"<getNumberOfCells()<<" cells and "<<_mesh[inew]->getNumberOfNodes()<<" nodes"<decrRef(); + } + std::multimap, pair > nodeMapping; + createNodeMapping(initial_collection, nodeMapping); + castMeshes(initial_collection.getFaceMesh(), this->getFaceMesh(),initial_collection, nodeMapping); + int nbOfGroups = initial_collection.getGroupMeshes(0).size(); + _groupMesh.resize(nbOfGroups); + for (int igroup=0; igroupgetGroupMeshes())[igroup], + initial_collection, nodeMapping); + for (int inew=0; inew<_topology->nbDomain(); inew++) + _groupMesh[igroup][inew]->setName((initial_collection.getGroupMeshes())[igroup][0]->getName()); + } +} + +void MESHCollection::createNodeMapping( MESHCollection& initialCollection, std::multimap,pair >& nodeMapping) +{ + for (int iold=0; ioldnbDomain();iold++) + { + std::map >, int > nodeClassifier; + for (int inode=0; inodegetNumberOfNodes(); inode++) + { + ParaMEDMEM::DataArrayDouble* coords = initialCollection.getMesh(iold)->getCoords(); + double* coordsPtr=coords->getPointer()+initialCollection.getMesh(iold)->getMeshDimension()*inode; + pair > tripleDouble =std::make_pair(coordsPtr[0],std::make_pair(coordsPtr[1],coordsPtr[2])); + nodeClassifier.insert(make_pair(tripleDouble,inode)); + + } + for (int inew=0; inew<_topology->nbDomain(); inew++) + { + for (int inode=0; inode<_mesh[inew]->getNumberOfNodes();inode++) + { + ParaMEDMEM::DataArrayDouble* coords = getMesh(inew)->getCoords(); + double* coordsPtr=coords->getPointer()+inode*getMesh(inew)->getMeshDimension(); + std::map >, int >::iterator iter= + nodeClassifier.find(make_pair(coordsPtr[0],make_pair(coordsPtr[1],coordsPtr[2]))); + nodeMapping.insert(make_pair(make_pair(iold,iter->second),make_pair(inew,inode))); + // cout <<"("<second<<")-->("<& meshesCastFrom,std::vector& meshesCastTo, MESHCollection& initialCollection,const multimap,pair >& nodeMapping) +{ + + //splitMeshes structure will contain the partition of + //the old faces on the new ones + // splitMeshes[4][2] contains the faces from old domain 2 + // that have to be added to domain 4 + + vector< vector > splitMeshes; + int newSize=_topology->nbDomain(); + splitMeshes.resize ( newSize ); + for (int inew=0;inew > newFaceMeshesIds(newSize); + for (int ielem=0;ielemgetNumberOfCells();ielem++) + { + std::vector nodes; + meshesCastFrom[iold]->getNodeIdsOfCell(ielem,nodes); + + map faces; + // cout<<"----------------"<,pair >::const_iterator MI; + pair myRange = nodeMapping.equal_range(make_pair(iold,nodes[inode])); + // cout << iold <<" " <second.first)==faces.end()) + faces[iter->second.first]=1; + else + faces[iter->second.first]++; + // cout<<"idomain" << iter->second.first<<" facemapping "<second.first]<::iterator iter=faces.begin(); + iter!=faces.end(); + iter++) + { + if (iter->second==nodes.size()) + newFaceMeshesIds[iter->first].push_back(ielem); + } + } + + //creating the splitMeshes from the face ids + for (int inew=0;inew<_topology->nbDomain();inew++) + { + // cout<<"nb faces "<buildPartOfMySelf(&newFaceMeshesIds[inew][0],&newFaceMeshesIds[inew][0]+newFaceMeshesIds[inew].size(),true)); + } + } + +// send/receive stuff + + +//recollecting the bits of splitMeshes to fuse them into one + meshesCastTo.resize(newSize); + for (int inew=0; inew < newSize;inew++) + { + vector myMeshes(meshesCastFrom.size()); + for (int iold=0; iold < meshesCastFrom.size();iold++) + { + myMeshes[iold]=splitMeshes[inew][iold]; + cout<<"number of nodes"<getNumberOfNodes()<getNumberOfCells()<getMeshDimension()<zipCoords(); + for (int iold=0; iold < meshesCastFrom.size();iold++) + splitMeshes[inew][iold]->decrRef(); + cout<<"creating face mesh #"<getNumberOfCells()<<" cells and "<getNumberOfNodes()<<" nodes"<read (filenamechar); + _driver_type = MedXML; + + } + catch(MEDMEM::MEDEXCEPTION){ + delete _driver; + try + { + _driver=new MESHCollectionMedAsciiDriver(this); + _driver->read (filenamechar); + _driver_type=MedAscii; + } + catch(MEDMEM::MEDEXCEPTION&) + { + delete _driver; + throw MEDMEM::MEDEXCEPTION("file does not comply with any recognized format"); + } + } + for ( int idomain = 0; idomain < _mesh.size(); ++idomain ) + if ( _mesh[idomain] && _mesh[idomain]->getNumberOfNodes() > 0 ) + _i_non_empty_mesh = idomain; +} + +/*! Constructing the MESH collection from selected domains of a distributed file + * + * \param filename - name of the master file containing the list of all the MED files + * \param domainSelector - selector of domains to load + */ +MESHCollection::MESHCollection(const string& filename, ParaDomainSelector& domainSelector) + : _topology(0), + _owns_topology(true), + _driver(0), + _domain_selector( domainSelector.nbProcs() > 1 ? & domainSelector : 0 ), + _i_non_empty_mesh(-1), + _driver_type(MEDPARTITIONER::Undefined), + _subdomain_boundary_creates(false), + _family_splitting(false), + _create_empty_groups(false) +{ + try + { + _driver=new MESHCollectionMedXMLDriver(this); + _driver->read ( (char*)filename.c_str(), _domain_selector ); + _driver_type = MedXML; + } + catch(MEDMEM::MEDEXCEPTION&) + { + delete _driver; + try + { + _driver=new MESHCollectionMedAsciiDriver(this); + _driver->read ( (char*)filename.c_str(), _domain_selector ); + _driver_type=MedAscii; + } + catch(MEDMEM::MEDEXCEPTION&) + { + delete _driver; + throw MEDMEM::MEDEXCEPTION("file does not comply with any recognized format"); + } + } + if ( isParallelMode() ) + // to know nb of cells on each proc to compute global cell ids from locally global + _domain_selector->gatherNbOf( getMesh() ); + + // find non-empty domain mesh + for ( int idomain = 0; idomain < _mesh.size(); ++idomain ) + if ( _mesh[idomain] && _mesh[idomain]->getNumberOfNodes() > 0 ) + _i_non_empty_mesh = idomain; +} + +/*! constructing the MESH collection from a sequential MED-file + * + * \param filename MED file + * \param meshname name of the mesh that is to be read + */ +MESHCollection::MESHCollection(const string& filename, const string& meshname) + : _name(meshname), + _topology(0), + _owns_topology(true), + _driver(0), + _domain_selector( 0 ), + _i_non_empty_mesh(-1), + _driver_type(MEDPARTITIONER::MedXML), + _subdomain_boundary_creates(false), + _family_splitting(false), + _create_empty_groups(false) +{ + char filenamechar[256]; + char meshnamechar[256]; + strcpy(filenamechar,filename.c_str()); + strcpy(meshnamechar,meshname.c_str()); + try // avoid memory leak in case of inexistent filename + { + retrieveDriver()->readSeq (filenamechar,meshnamechar); + } + catch ( MED_EXCEPTION& e ) + { + if ( _driver ) delete _driver; _driver=0; + throw e; + } + if ( _mesh[0] && _mesh[0]->getNumberOfNodes() > 0 ) + _i_non_empty_mesh = 0; +} + +MESHCollection::~MESHCollection() +{ + for (int i=0; i<_mesh.size();i++) + if (_mesh[i]!=0) {/*delete*/ _mesh[i]->decrRef(); } + for (int i=0; i<_faceMesh.size();i++) + if (_mesh[i]!=0) {/*delete*/ _faceMesh[i]->decrRef(); } + for (int i=0; i<_groupMesh.size();i++) + for (int j=0; j<_groupMesh[i].size(); j++) + _groupMesh[i][j]->decrRef(); + for (int i=0; i<_connect_zones.size();i++) + if (_connect_zones[i]!=0) {delete _connect_zones[i];} + if (_driver !=0) {delete _driver; _driver=0;} + if (_topology!=0 && _owns_topology) {delete _topology; _topology=0;} +} + + +/*! constructing the MESH collection from a file + * + * The method creates as many MED-files as there are domains in the + * collection. It also creates a master file that lists all the MED files. + * The MED files created in ths manner contain joints that describe the + * connectivity between subdomains. + * + * \param filename name of the master file that will contain the list of the MED files + * + */ +void MESHCollection::write(const string& filename) +{ + //building the connect zones necessary for writing joints + cout<<"Building Connect Zones"<nbDomain()>1) +// buildConnectZones(); + cout <<"End of connect zones building"<write (filenamechar, _domain_selector); +} + +/*! creates or gets the link to the collection driver + */ +MESHCollectionDriver* MESHCollection::retrieveDriver() +{ + if (_driver==0) + { + switch(_driver_type) + { + case MedXML: + _driver=new MESHCollectionMedXMLDriver(this); + break; + case MedAscii: + _driver=new MESHCollectionMedAsciiDriver(this); + break; + default: + throw MEDMEM::MEDEXCEPTION("Unrecognized driver"); + } + } + + return _driver; +} + + +/*! gets an existing driver + * + */ +MESHCollectionDriver* MESHCollection::getDriver() const +{ + return _driver; +} +// /*! retrieves the mesh dimension*/ + int MESHCollection::getMeshDimension() const + { + return _i_non_empty_mesh < 0 ? -1 : _mesh[_i_non_empty_mesh]->getMeshDimension(); + } + + vector& MESHCollection::getMesh() + { + return _mesh; + } + vector& MESHCollection::getFaceMesh() + { + return _faceMesh; + } +vector >& MESHCollection::getGroupMeshes() + { + return _groupMesh; + } +ParaMEDMEM::MEDCouplingUMesh* MESHCollection::getMesh(int idomain) +{ + return _mesh[idomain]; +} + +ParaMEDMEM::MEDCouplingUMesh* MESHCollection::getFaceMesh(int idomain) +{ + return _faceMesh[idomain]; +} +vector& MESHCollection::getGroupMeshes(int idomain) +{ + return _groupMesh[idomain]; +} + vector& MESHCollection::getCZ() + { + return _connect_zones; + } + +Topology* MESHCollection::getTopology() const + { + return _topology; + } + +void MESHCollection::setTopology(Topology* topo) + { + if (_topology!=0) + { + throw MEDMEM::MEDEXCEPTION("Erreur : topology is already set"); + } + else + _topology = topo; + } + + +/*! Method creating the cell graph + * + * \param array returns the pointer to the structure that contains the graph + * \param edgeweight returns the pointer to the table that contains the edgeweights + * (only used if indivisible regions are required) + */ + +void MESHCollection::buildCellGraph(MEDPARTITIONER::MEDSKYLINEARRAY* & array,int *& edgeweights ) +{ + + multimap< int, int > node2cell; + multimap< int, int > cell2cell; + + //looking for reverse nodal connectivity i global numbering + for (int idomain=0; idomain<_topology->nbDomain();idomain++) + { + ParaMEDMEM::DataArrayInt* index=ParaMEDMEM::DataArrayInt::New(); + ParaMEDMEM::DataArrayInt* revConn=ParaMEDMEM::DataArrayInt::New(); + _mesh[idomain]->getReverseNodalConnectivity(revConn,index); + int* index_ptr=index->getPointer(); + ParaMEDMEM::DataArrayInt* globalRevConn=revConn->deepCpy(); + _topology->convertCellToGlobal(idomain, + revConn->getPointer(), + index_ptr[_mesh[idomain]->getNumberOfNodes()], + globalRevConn->getPointer()); + + int min=10000000,max=-1; + for (int i=0; i< index_ptr[_mesh[idomain]->getNumberOfNodes()];i++) + { + if (revConn->getPointer()[i]getPointer()[i]; + if (revConn->getPointer()[i]>max) max=revConn->getPointer()[i]; + } + cout <<"min:"<getNumberOfNodes()]; + _topology->getNodeList(idomain,globalNodeIds); + + int* globalRevConnPtr=globalRevConn->getPointer(); + for (int i=0; i<_mesh[idomain]->getNumberOfNodes();i++) + { + for (int icell=index_ptr[i]; icelldecrRef(); + revConn->decrRef(); + index->decrRef(); + delete[]globalNodeIds; + } + + //creating graph arcs (cell to cell relations) + //arcs are stored in terms of (index,value) notation + // 0 3 5 6 6 + // 1 2 3 2 3 3 + // means 6 arcs (0,1), (0,2), (0,3), (1,2), (1,3), (2,3) + // in present version arcs are not doubled but reflexive (1,1) arcs are present for each cell + + for (int inode=0; inode<_topology->nbNodes();inode++) + { + typedef multimap::const_iterator MI; + pair myRange = node2cell.equal_range(inode); + for (MI cell1=myRange.first; cell1!=myRange.second; cell1++) + { + for (MI cell2 = myRange.first; cell2!=myRange.second; cell2++) + { + if (cell1->second!=cell2->second) cell2cell.insert(make_pair(cell1->second,cell2->second)); + } + } + } + + //filling up index and value to create skylinearray structure + + multimap::const_iterator iter; + + iter=cell2cell.begin(); + vector index,value; + index.push_back(0); + while (iter != cell2cell.end()) + { + multimap::const_iterator next_iter = cell2cell.upper_bound(iter->first); + int size=0; + for (multimap::const_iterator current_iter=iter; current_iter!=next_iter; current_iter++) + { + value.push_back(current_iter->second); + size++; + } + index.push_back(index[index.size()-1]+size); + + + iter=next_iter; + } + + array=new MEDPARTITIONER::MEDSKYLINEARRAY(index,value); + cout<<"taille du graphe créé "<getNumberOf()<nbDomain(); i++) +// { +// cell_number+=_topology->getCellNumber(i); +// node_number+=_topology->getNodeNumber(i); +// } +// //list of cells for a given node +// //vector< vector > node2cell(node_number); +// map< int, vector > node2cell; + +// //list of nodes for a given cell +// //vector< vector > cell2node(cell_number); +// map< int, vector > cell2node; + +// // map type_cell_list; + +// //tagging for the indivisible regions +// int* indivisible_tag=0; +// bool has_indivisible_regions=false; +// // if (!_indivisible_regions.empty()) +// // { +// // has_indivisible_regions=true; +// // indivisible_tag=new int[_topology->nbCells()]; +// // treatIndivisibleRegions(indivisible_tag); +// // } + +// fillGlobalConnectivity(node2cell, cell2node ); + +// cout << "beginning of skyline creation"<nbCells()]; +// int** temp=new int*[_topology->nbCells()]; +// int** temp_edgeweight=0; +// // if (has_indivisible_regions) +// // temp_edgeweight=new int*[_topology->nbCells()]; + +// int cell_glob_shift = 0; + +// // Get connection to cells on other procs +// multimap< int, int > loc2dist; // global cell ids on this proc -> other proc cells +// if ( isParallelMode() ) +// { +// cell_glob_shift = _domain_selector->getProcShift(); + +// set loc_domains; // domains on this proc +// for ( int idom = 0; idom < _mesh.size(); ++idom ) +// if ( _mesh[ idom ] ) +// loc_domains.insert( idom ); + +// for ( int idom = 0; idom < _mesh.size(); ++idom ) +// { +// if ( !_mesh[idom] ) continue; +// vector loc2glob_corr; // pairs of corresponding cells (loc_loc & glob_dist) +// retrieveDriver()->readLoc2GlobCellConnect(idom, loc_domains, _domain_selector, loc2glob_corr); +// //MEDMEM::STRING out; +// for ( int i = 0; i < loc2glob_corr.size(); i += 2 ) +// { +// int glob_here = _topology->convertCellToGlobal(idom,loc2glob_corr[i]); +// int glob_there = loc2glob_corr[i+1]; +// loc2dist.insert ( make_pair( glob_here, glob_there)); +// //out << glob_here << "-" << glob_there << " "; +// } +// //cout << "\nRank " << _domain_selector->rank() << ": BndCZ: " << out << endl; +// } +// } + +// //going across all cells + +// map cells_neighbours; +// for (int i=0; i< _topology->nbCells(); i++) +// { + + +// vector cells(50); + +// // /*// cout << "size cell2node "<::const_iterator iternode=cell2node[i+1].begin(); +// // iternode!=cell2node[i+1].end(); +// // iternode++) +// // { +// // int nodeid=*iternode; +// // // cout << "size node2cell "<::const_iterator iter=node2cell[nodeid].begin(); +// // iter != node2cell[nodeid].end(); +// // iter++) +// // cells_neighbours[*iter]++; +// // }*/ + +// for (int inode=0; inode< cell2node[i+1].size(); inode++) +// { +// int nodeid=cell2node[i+1][inode]; +// // cout << "size node2cell "<::const_iterator iter=cells_neighbours.begin(); iter != cells_neighbours.end(); iter++) +// { +// if (iter->second >= dimension && iter->first != i+1) +// { +// cells.push_back(iter->first + cell_glob_shift); +// // cells[isize++]=iter->first; +// } +// } +// // add neighbour cells from distant domains +// multimap< int, int >::iterator loc_dist = loc2dist.find( i+1 ); +// for (; loc_dist!=loc2dist.end() && loc_dist->first==( i+1 ); ++loc_dist ) +// cells.push_back( loc_dist->second ); + +// size[i]=cells.size(); +// // size[i]=isize; + +// //cout << cells.size()<::const_iterator iter=cells.begin(); iter!=cells.end();iter++) +// //for(int j=0; jfirst-1]; +// // int tag2 = indivisible_tag[*iter-1]; +// // if (tag1==tag2 && tag1!=0) +// // temp_edgeweight[i][itemp]=_topology->nbCells()*100000; +// // else +// // temp_edgeweight[i][itemp]=1; +// // } +// itemp++; +// } +// cells_neighbours.clear(); +// } +// cout <<"end of graph definition"<nbCells()+1]; +// index[0]=1; +// for (int i=0; i<_topology->nbCells(); i++) +// index[i+1]=index[i]+size[i]; + +// node2cell.clear(); +// cell2node.clear(); +// // if (indivisible_tag!=0) delete [] indivisible_tag; + +// //SKYLINEARRAY structure holding the cell graph +// array= new MEDMEM::MEDSKYLINEARRAY(_topology->nbCells(),index[_topology->nbCells()]-index[0]); +// array->setIndex(index); + +// for (int i=0; i<_topology->nbCells(); i++) +// { +// array->setI(i+1,temp[i]); +// delete[]temp[i]; +// } + +// {// DEBUG +// // MEDMEM::STRING out; +// // const int* index= array->getIndex(); +// // const int* value =array->getValue(); +// // out << "\nRank " << (_domain_selector?_domain_selector->rank():0) << ": Index: "; +// // for ( int i = 0; i <= array->getNumberOf(); ++i ) +// // out << index[i] << " "; +// // out << "\nRank " << (_domain_selector?_domain_selector->rank():0) << ": Value: "; +// // for ( int i = 0; i < array->getLength(); ++i ) +// // out << value[i] << " "; +// // cout << out << endl; +// } +// // if (has_indivisible_regions) +// // { +// // edgeweights=new int[array->getLength()]; +// // for (int i=0; i<_topology->nbCells(); i++) +// // { +// // for (int j=index[i]; j0"); + MEDPARTITIONER::MEDSKYLINEARRAY* array=0; + int* edgeweights=0; + + cout<<"Building cell graph"; + // if ( _domain_selector ) + // buildCellGraphParallel(array,edgeweights); + // else + buildCellGraph(array,edgeweights); + Graph* cellGraph; + switch (split) + { + case Graph::METIS: +#ifdef ENABLE_METIS + cellGraph=(Graph*)(new METISGraph(array,edgeweights)); +#else + throw MEDMEM::MEDEXCEPTION("METIS Graph is not available. Check your products, please."); +#endif + break; + case Graph::SCOTCH: +#ifdef ENABLE_SCOTCH + cellGraph=(Graph*)(new SCOTCHGraph(array,edgeweights)); +#else + throw MEDMEM::MEDEXCEPTION("SCOTCH Graph is not available. Check your products, please."); +#endif + break; + } + + //!user-defined weights + if (user_edge_weights!=0) + cellGraph->setEdgesWeights(user_edge_weights); + if (user_vertices_weights!=0) + cellGraph->setVerticesWeights(user_vertices_weights); + + cout<<"Partitioning graph"; + cellGraph->partGraph(nbdomain,options_string,_domain_selector); + + // DEBUG + // MEDMEM::STRING out("RESULT GRAPH #"); + // out << (_domain_selector?_domain_selector->rank():0) << ": "; + // const int* part = cellGraph->getPart(); + // int n = cellGraph->nbVertices(); + // for ( int e=0; e < n; ++e ) + // out << part[e] <<" "; + // cout << out << endl; + + + cout<<"Building new topology"; + //cellGraph is a shared pointer + Topology* topology = new ParallelTopology (cellGraph, nbdomain, getMeshDimension()); + + //cleaning + if (edgeweights!=0) delete[] edgeweights; + // if (array!=0) delete array; + cout<<"End of partition creation"; + delete cellGraph; + + return topology; +} + +/*! Creates a topology for a partition specified by the user + * + * \param table user-specified partition (for each cell contains the domain number from 0 to n-1) + * + * returns a topology based on the new partition + */ +Topology* MESHCollection::createPartition(const int* partition) +{ + MEDPARTITIONER::MEDSKYLINEARRAY* array=0; + int* edgeweights=0; + + buildCellGraph(array,edgeweights); + Graph* cellGraph; + set domains; + for (int i=0; i<_topology->nbCells(); i++) + { + domains.insert(partition[i]); + } + int nbdomain=domains.size(); + + cellGraph=(Graph*)(new UserGraph(array, partition, _topology->nbCells())); + + //cellGraph is a shared pointer + Topology* topology = new ParallelTopology (cellGraph, nbdomain, getMeshDimension()); + + // if (array!=0) delete array; + delete cellGraph; + return topology; +} + + +/*! building Connect Zones for storing the informations + * of the connectivity + * + * The connect zones are created for every domain that has common nodes with + * domain \a idomain + * + * \param idomain domain number for which the connect zones are created + * */ + +// void MESHCollection::buildConnectZones(int idomain) +// { +// // constructing node/node correspondencies +// vector node_node_correspondency; +// node_node_correspondency.resize(_topology->nbDomain()); + +// cout << "Computing node/node corresp"<computeNodeNodeCorrespondencies(idomain, node_node_correspondency ); + +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// { +// // on regarde si une correspondance noeud/noeud a été trouvée +// // entre idomain et idistant +// // si oui, on crée une connectzone +// if (node_node_correspondency[idistant]!=0) +// { +// MEDMEM::CONNECTZONE* cz= new MEDMEM::CONNECTZONE(); +// cz->setLocalMesh(_mesh[idomain]); +// cz->setDistantMesh(_mesh[idistant]); +// cz->setLocalDomainNumber(idomain); +// cz->setDistantDomainNumber(idistant); +// cz-> setName ("Connect zone defined by PARTITIONER"); +// cz->setNodeCorresp(node_node_correspondency[idistant]); +// _connect_zones.push_back(cz); +// } +// } +// cout << "Computing node/node corresp"< cell_cell_correspondency; +// cell_cell_correspondency.resize(_topology->nbDomain()); +// _topology->computeCellCellCorrespondencies(idomain, cell_cell_correspondency, _cell_graph.get()); + +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// { +// //the connect zone has been created by the node/node computation +// if (cell_cell_correspondency[idistant]!=0) +// { +// MEDMEM::CONNECTZONE* cz=0; +// for (int icz=0; icz<_connect_zones.size();icz++) +// if (_connect_zones[icz]->getLocalDomainNumber()==idomain && +// _connect_zones[icz]->getDistantDomainNumber()==idistant) +// cz = _connect_zones[icz]; +// if (cz!=0) +// cz->setEntityCorresp(MED_EN::MED_CELL,MED_EN::MED_CELL, cell_cell_correspondency[idistant]); +// else +// throw MEDEXCEPTION("MESHCollection::buildConnectZones() -A connect zone should exist"); +// //delete cell_cell_correspondency[idistant]; +// } + +// } +// } + +// //================================================================================ +// /*! +// * \brief Adds a group of joint faces +// * \param loc_face_ids - local numbers of faces +// * \param idomian - domain index where faces are local +// * \param idistant - the other domain index +// */ +// //================================================================================ + +// void MESHCollection::addJointGroup(const std::vector& loc_face_ids, int idomain, int idistant) +// { +// MEDMEM::MESHING* meshing = dynamic_cast (_mesh[idomain]); +// MED_EN::medEntityMesh constituent_entity = getSubEntity(); + +// MEDMEM::STRING jointname("joint_"); +// jointname< groups = meshing->getGroups( constituent_entity ); +// // for ( int g = 0; g < groups.size(); ++g ) +// // if ( groups[g]->getName() == jointname.str() ) +// // { +// // joint_group = groups[g]; +// // break; +// // } +// // assure uniqueness of group name +// bool unique = false; +// vector groups = meshing->getGroups( constituent_entity ); +// do +// { +// unique = true; +// for ( int g = 0; unique && g < groups.size(); ++g ) +// unique = ( groups[g]->getName() != jointname ); +// if ( !unique ) +// jointname << "_" << idomain+1; +// } +// while ( !unique ); +// joint_group->setName(jointname); +// joint_group->setMesh(meshing); +// joint_group->setEntity(constituent_entity); +// map > joint_types; + +// int nbfaces = loc_face_ids.size(); +// for (int i=0; igetElementType(constituent_entity,loc_face_ids[i]); +// joint_types[type].push_back(loc_face_ids[i]); +// } +// joint_group->setNumberOfGeometricType(joint_types.size()); +// MED_EN::medGeometryElement* types=new MED_EN::medGeometryElement[joint_types.size()]; +// int* nb_in_types=new int[joint_types.size()]; +// int* group_index=new int[joint_types.size()+1]; + +// group_index[0]=1; +// int itype=0; +// int iface =0; +// int* group_value=new int[nbfaces]; +// for (map >::const_iterator iterj=joint_types.begin(); +// iterj != joint_types.end(); +// iterj++) +// { +// nb_in_types[itype]=(iterj->second).size(); +// types[itype]=iterj->first; +// itype++; +// group_index[itype]=group_index[itype-1]+(iterj->second).size(); +// for (int i=0; i< (iterj->second).size(); i++) +// group_value[iface++]=(iterj->second)[i]; +// } +// joint_group->setGeometricType(types); +// joint_group->setNumberOfElements(nb_in_types); +// joint_group->setNumber(group_index, group_value, /*shallowCopy=*/true); +// delete[] types; +// delete[] nb_in_types; + +// if ( joint_group == tmp_grp ) +// meshing->addGroup(*tmp_grp); +// tmp_grp->removeReference(); +// } + +// /*! building Connect Zones for storing the informations +// * of the connectivity +// * */ + +// void MESHCollection::buildConnectZones() +// { +// vector > > face_map(_topology->nbDomain()); +// map< pair, MEDMEM::MEDSKYLINEARRAY*> cell_corresp_here; + +// MED_EN::medEntityMesh constituent_entity = getSubEntity(); + +// if ( isParallelMode() ) +// { +// buildConnectZonesBetweenProcs(face_map, cell_corresp_here); +// } + +// cout << "Computing node/node corresp"<nbDomain(); idomain++) +// { + +// // constructing node/node correspondencies +// vector node_node_correspondency(_topology->nbDomain()); +// _topology->computeNodeNodeCorrespondencies(idomain, node_node_correspondency ); + +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// { +// // on regarde si une correspondance noeud/noeud a été trouvée +// // entre idomain et idistant +// // si oui, on crée une connectzone +// if (node_node_correspondency[idistant]!=0) +// { +// MEDMEM::CONNECTZONE* cz= new MEDMEM::CONNECTZONE(); +// cz->setLocalMesh(_mesh[idomain]); +// cz->setDistantMesh(_mesh[idistant]); +// cz->setLocalDomainNumber(idomain); +// cz->setDistantDomainNumber(idistant); +// cz-> setName ("Connect zone defined by PARTITIONER"); +// cz->setNodeCorresp(node_node_correspondency[idistant]); +// _connect_zones.push_back(cz); +// } +// } +// } +// cout << "Computing face corresp"<getFaceNumber()+1; +// //int global_face_id = _topology->getMaxGlobalFace()+1; + +// map , vector > faces_in_joint; + +// if ( !isParallelMode() ) +// // taking faces that are already present in the mesh into account +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// { +// getFaces(idomain,face_map[idomain]); +// } + +// // creating faces that are located at the interface between +// // subdomains + +// vector nb_added_groups( _topology->nbDomain(), 0 ); + +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// { +// vector cell_cell_correspondency( _topology->nbDomain() ); +// if ( !isParallelMode() ) +// _topology->computeCellCellCorrespondencies(idomain, cell_cell_correspondency, _cell_graph.get()); + +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// { +// if (idistant <= idomain) continue; + +// MEDMEM::MEDSKYLINEARRAY* cell_correspondency = 0; +// if ( isParallelMode() ) +// cell_correspondency = cell_corresp_here[ make_pair (idomain,idistant)]; +// else +// cell_correspondency = cell_cell_correspondency[idistant]; + +// //the connect zone has been created by the node/node computation + +// if ( cell_correspondency ) +// { +// int nbcells = cell_correspondency->getNumberOf(); +// const int* index = cell_correspondency->getIndex(); +// const int* value = cell_correspondency->getValue(); +// if ( isParallelMode() ) +// global_face_id = _domain_selector->getFisrtGlobalIdOfSubentity( idomain, idistant ); + +// for (int ilocal=0; ilocalgetType()].push_back(face); +// MEDPARTITIONER_FaceModel* face2 = getCommonFace(idistant,distant_id,idomain, ilocal+1,global_face_id); +// face_map[idistant][face->getType()].push_back(face2); +// faces_in_joint[make_pair(idomain,idistant)].push_back(global_face_id); +// global_face_id++; +// } +// } +// } + +// } +// //cleaning up +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// delete cell_cell_correspondency[idistant]; +// } + + +// _topology->recreateFaceMapping(face_map); + +// //transforming the face_map into a constituent entity connectivity +// for (int idomain=0; idomain< _topology->nbDomain();idomain++) +// { +// int nbtypes = face_map[idomain].size(); +// vector types; +// vector nb_elems; +// vector conn; + +// MEDMEM::MESHING* meshing = dynamic_cast (_mesh[idomain]); +// if ( !meshing->getConnectivityptr() ) +// continue; // no cells in idomain + +// for (map >::const_iterator iter= face_map[idomain].begin(); +// iter != face_map[idomain].end(); iter ++) +// { +// types.push_back(iter->first); +// int nb_elem_in_type = (iter->second).size(); +// nb_elems.push_back(nb_elem_in_type); +// int nb_node_per_type=(iter->first)%100; +// int* connectivity= new int [nb_node_per_type*nb_elem_in_type]; +// for (int ielem=0; ielemsecond)[ielem])[inode]; +// } +// conn.push_back(connectivity); + +// } +// //setting the faces in the mesh +// meshing->setNumberOfTypes(nbtypes,constituent_entity); +// meshing->setTypes(&types[0],constituent_entity); +// meshing->setNumberOfElements(&nb_elems[0],constituent_entity); + +// for (int itype=0; itypesetConnectivity(conn[itype], constituent_entity, types[itype]); +// delete[]conn[itype]; +// } +// for (int idistant =0; idistant<_topology->nbDomain(); idistant++) +// { +// map , vector >::iterator iter; +// iter = faces_in_joint.find(make_pair(idomain,idistant)); +// if (iter == faces_in_joint.end()) +// { +// iter = faces_in_joint.find (make_pair(idistant,idomain)); +// if (iter == faces_in_joint.end()) +// continue; +// } + +// int nbfaces = (iter->second).size(); +// vector face_joint(nbfaces*2); +// MEDMEM::CONNECTZONE* cz=0; +// for (int icz=0; icz<_connect_zones.size();icz++) +// if (_connect_zones[icz]->getLocalDomainNumber()==idomain && +// _connect_zones[icz]->getDistantDomainNumber()==idistant) +// cz = _connect_zones[icz]; + +// int nbtotalfaces= _topology->getFaceNumber(idomain); + +// //creating arrays for the MEDSKYLINEARRAY structure containing the joint +// int* index =new int[nbtotalfaces+1]; +// for (int i=0; i faces; +// vector local_faces( nbfaces ); +// for (int iface=0; ifacesecond)[iface]; +// int localid=_topology->convertGlobalFace(iglobal,idomain); +// int distantid=_topology->convertGlobalFace(iglobal,idistant); +// faces.insert(make_pair(localid,distantid)); +// local_faces[iface]=localid; +// } + +// int iloc=0; +// index[0]=1; +// for (map::const_iterator iter=faces.begin(); +// iter != faces.end(); +// iter++) +// { +// index[iter->first]=1; +// value[iloc++]=iter->second; +// } + +// for (int i=0; isetEntityCorresp(constituent_entity,constituent_entity,skarray); +// else +// throw MEDEXCEPTION("MESHCollection::buildConnectZones() -A connect zone should exist"); +// // Creating a group of the faces constituting the joint +// addJointGroup( local_faces, idomain, idistant ); +// nb_added_groups[ idomain ]++; +// } +// } + +// if ( isParallelMode() ) +// { +// // Now all faces have got local ids and we can receive local ids from other procs. +// // Set face/face data to zones with other procs and create a group +// for (int icz=0; icz<_connect_zones.size();icz++) +// { +// MEDMEM::CONNECTZONE* cz=_connect_zones[icz]; +// if ( _domain_selector->isMyDomain( cz->getDistantDomainNumber()) ) continue; + +// int glob_id = _domain_selector->getFisrtGlobalIdOfSubentity( cz->getLocalDomainNumber(), +// cz->getDistantDomainNumber()); +// int nb_cz_faces = _domain_selector->getNbCellPairs( cz->getDistantDomainNumber(), +// cz->getLocalDomainNumber()); +// vector< int > loc_ids_here( nb_cz_faces ); +// for ( int i = 0; i < nb_cz_faces; ++i ) +// loc_ids_here[i] = _topology->convertGlobalFace(glob_id++,cz->getLocalDomainNumber()); + +// int* loc_ids_dist = _domain_selector->exchangeSubentityIds( cz->getLocalDomainNumber(), +// cz->getDistantDomainNumber(), +// loc_ids_here ); +// int nb_faces_here= _topology->getFaceNumber(cz->getLocalDomainNumber()); +// int* face_index = new int[ nb_faces_here+1 ]; +// face_index[0]=1; +// for ( int loc_id = 0, i = 0; loc_id < nb_faces_here; ++loc_id) +// { +// face_index[ loc_id+1 ] = face_index[ loc_id ]; +// if ( i < loc_ids_here.size() && loc_ids_here[i] == loc_id+1 ) +// { +// face_index[ loc_id+1 ]++; +// i++; +// } +// } +// MEDMEM::MEDSKYLINEARRAY* skarray= +// new MEDMEM::MEDSKYLINEARRAY(nb_faces_here, nb_cz_faces, face_index, loc_ids_dist, true); +// cz->setEntityCorresp(constituent_entity,constituent_entity,skarray); + +// addJointGroup( loc_ids_here, cz->getLocalDomainNumber(), cz->getDistantDomainNumber()); +// nb_added_groups[ cz->getLocalDomainNumber() ]++; +// } +// } + +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// { +// // delete face_map +// for (map >::const_iterator iter= face_map[idomain].begin(); +// iter != face_map[idomain].end(); iter ++) +// for (int i=0; i<(iter->second).size();i++) +// delete (iter->second)[i]; + +// if ( nb_added_groups[ idomain ] > 0 && +// _mesh[idomain]->getNumberOfFamilies( constituent_entity ) > 0 ) +// // needed because if there were face families before, driver won't +// // create families from just added groups (see MEDMEM_MedMeshDriver22.cxx:3330), +// // actually it is a bug of driver - it must check presence of groups in families +// _mesh[idomain]->createFamilies(); +// } +// } + +// if ( isParallelMode() ) +// // Excange info on types of constituent_entity needed while writing joints +// // to get ids local in geom type for distant procs +// _domain_selector->gatherEntityTypesInfo( _mesh, constituent_entity ); + +// cout << "Computing cell/cell corresp"<nbDomain();idomain++) +// { +// vector cell_cell_correspondency( _topology->nbDomain() ); +// if ( !isParallelMode() ) +// _topology->computeCellCellCorrespondencies(idomain,cell_cell_correspondency,_cell_graph.get()); + +// for (int idistant=0; idistant< _topology->nbDomain(); idistant++) +// { +// MEDMEM::MEDSKYLINEARRAY* cell_correspondency = 0; +// if ( isParallelMode() ) +// cell_correspondency = cell_corresp_here[ make_pair (idomain,idistant)]; +// else +// cell_correspondency = cell_cell_correspondency[idistant]; + +// //the connect zone has been created by the node/node computation +// if ( cell_correspondency ) +// { +// MEDMEM::CONNECTZONE* cz=0; +// for (int icz=0; icz<_connect_zones.size();icz++) +// if (_connect_zones[icz]->getLocalDomainNumber()==idomain && +// _connect_zones[icz]->getDistantDomainNumber()==idistant) +// cz = _connect_zones[icz]; +// if (cz!=0) +// cz->setEntityCorresp(MED_EN::MED_CELL,MED_EN::MED_CELL, cell_correspondency); +// else +// throw MEDEXCEPTION("MESHCollection::buildConnectZones() -A connect zone should exist"); +// } +// } +// } +// } + +// /*! building Connect Zones for storing the informations +// * of the connectivity in the parallel mode +// * */ + +// void MESHCollection::buildConnectZonesBetweenProcs(TGeom2FacesByDomian & face_map, +// map< pair, MEDMEM::MEDSKYLINEARRAY*> & cell_cell_correspondency_here) +// { +// using namespace MED_EN; + +// // graph over all procs +// auto_ptr global_graph( _domain_selector->gatherGraph( _cell_graph.get() )); + +// vector< vector< JointExchangeData > > joints_of_domain( _topology->nbDomain() ); + +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// { +// if ( !_domain_selector->isMyDomain( idomain )) continue; + +// vector< JointExchangeData > & joints = joints_of_domain[ idomain ]; +// joints.resize( _topology->nbDomain() ); + +// // Find corresponding cells on other procs + +// const int* gra_index = global_graph->getGraph()->getIndex(); +// const int* gra_value = global_graph->getGraph()->getValue(); +// const int* partition = global_graph->getPart(); +// const int dj = gra_index[0]; + +// vector< int > glob_cells_here( _topology->getCellNumber( idomain )); +// _topology->getCellList( idomain, & glob_cells_here[0]); +// for ( int loc_here = 0; loc_here < glob_cells_here.size(); ++loc_here ) +// { +// int glob_here = glob_cells_here[ loc_here ]; +// for ( int j = gra_index[ glob_here-1 ]; j < gra_index[ glob_here ]; ++j ) +// { +// int glob_neighbor = gra_value[ j-dj ]; +// int neighbor_dom = partition[ glob_neighbor-1 ]; +// if ( neighbor_dom == idomain ) continue; + +// if ( _domain_selector->isMyDomain( neighbor_dom )) +// { +// joints[ neighbor_dom ].addCellCorrespondence +// (_mesh[idomain], neighbor_dom, idomain, glob_neighbor, glob_here, loc_here + 1, +// _topology->convertGlobalCell(glob_neighbor).second ); +// } +// else +// { +// joints[ neighbor_dom ].addCellCorrespondence +// (_mesh[idomain], neighbor_dom, idomain, glob_neighbor, glob_here, loc_here + 1 ); +// } +// } +// } +// } +// global_graph.reset(); // free memory + +// // set joints in a queue to exchange +// typedef map< int, JointExchangeData* > TOrderedJoints; +// TOrderedJoints queue; +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// { +// if ( !_domain_selector->isMyDomain( idomain )) continue; + +// vector< JointExchangeData > & joints = joints_of_domain[ idomain ]; +// for (int idist=0; idist<_topology->nbDomain(); ++idist ) +// { +// JointExchangeData& joint = joints[idist]; + +// int nb_cell_pairs = joint.nbCellPairs(); +// if ( nb_cell_pairs == 0 ) +// continue; +// else +// _domain_selector->setNbCellPairs( nb_cell_pairs, idist, idomain ); + +// joint.setMeshes( idist, _mesh[idist], idomain, _mesh[idomain] ); + +// if ( _domain_selector->isMyDomain( idist )) +// { +// // a joint on this proc +// cell_cell_correspondency_here[ make_pair( idomain, idist )] = joint.makeCellCorrespArray(); +// } +// else +// { +// // a joint with distant proc +// joint.setConnectivity( & ((MEDMEM::MeshFuse*)_mesh[idomain])->getNodeNumbers()[0] ); +// int order = _domain_selector->jointId( idomain, idist ); +// queue[ order ] = & joint; +// } +// } +// } +// // gather info on cell geom types needed to exchange joints +// _domain_selector->gatherEntityTypesInfo( _mesh, MED_EN::MED_CELL ); + +// // gather info on nb of sub-entities to compute their global numbers for joints +// _domain_selector->gatherNbOf( getSubEntity(), _mesh ); +// _domain_selector->gatherNbCellPairs(); +// if ( _subdomain_boundary_creates ) +// { +// // taking faces that are already present in the mesh into account +// for (int idomain=0; idomain<_topology->nbDomain(); idomain++) +// if ( _domain_selector->isMyDomain( idomain )) +// getFaces(idomain,face_map[idomain]); +// } +// else +// { +// face_map.clear(); // mark for the joint not to create faces +// } + +// // exchange joint data with other procs and make CONNECTZONEs +// TOrderedJoints::iterator ord_joint = queue.begin(); +// for ( ; ord_joint != queue.end(); ++ord_joint ) +// { +// JointExchangeData* joint = ord_joint->second; + +// _domain_selector->exchangeJoint( joint ); +// if ( _subdomain_boundary_creates ) +// { +// int first_sub_id = _domain_selector->getFisrtGlobalIdOfSubentity( joint->localDomain(), +// joint->distantDomain() ); +// joint->setFisrtGlobalIdOfSubentity( first_sub_id ); +// } +// _connect_zones.push_back ( joint->makeConnectZone( face_map )); +// } +// } + +// /*! projects old collection families on new collection families +// */ +// void MESHCollection::castFamilies(const MESHCollection& old_collection) +// { +// vector > element_array (_topology->nbDomain()); + +// //loop on old domains to create groups out of the existing families +// if (_family_splitting) +// for (int idomain=0; idomain < old_collection._topology->nbDomain(); idomain++) +// old_collection.getMesh(idomain)->createGroups(); + +// //definition of the entities array which +// //defines the entities over which the information is cast +// MED_EN::medEntityMesh entities[3]; +// entities[0]=MED_EN::MED_NODE; +// entities[1]=getSubEntity(); +// entities[2]=MED_EN::MED_CELL; + +// for (int ientity=0; ientity<=2;ientity++) +// { + +// //int nbgroups = old_collection.getMesh(0)->getNumberOfGroups(entities[ientity]); + +// map > group_map; +// for (int idomain=0; idomain < old_collection._topology->nbDomain(); idomain++) +// { +// if ( !old_collection.getMesh(idomain) ) continue; +// for (int igroup=0; igroupgetNumberOfGroups(entities[ientity]); igroup++) +// { +// //retrieves a group +// MEDMEM::GROUP* group = (old_collection.getMesh(idomain)->getGroups(entities[ientity]))[igroup]; +// //increments the number of groups if it is a new group +// //if (group_map.find(group->getName())==group_map.end()) + +// group_map[group->getName()].insert(idomain); +// // group_map.insert(make_pair(group->getName(), idomain); + +// } +// } +// int nbgroups=group_map.size(); +// vector igroupold(old_collection._topology->nbDomain(),0); +// map >::const_iterator iter=group_map.begin(); + +// for (int igroup=0; igroup old_supports(old_collection._topology->nbDomain()); +// string group_name = iter->first; +// iter++; + +// //parameters stored for passing group description +// // from the old meshes to the new ones + +// for (int idomain=0; idomain < old_collection._topology->nbDomain(); idomain++) +// { +// // for (set::iterator iter=group_map[group_name].begin(); iter!=group_map[group_name].end(); iter++) +// // cout << *iter<<" "; +// // cout <getGroups(entities[ientity]))[igroupold[idomain]]; +// old_supports[idomain] = static_cast (group); +// igroupold[idomain]++; +// } + +// vector new_groups(_topology->nbDomain()); +// vector new_supports(_topology->nbDomain()); +// for (int i=0; i<_topology->nbDomain(); i++) +// { +// new_groups[i]=new MEDMEM::GROUP(); +// new_supports[i]=static_cast(new_groups[i]); +// } +// castSupport(old_collection,old_supports,new_supports); + +// //creating new groups from the previous list of elements +// for (int idomain=0; idomain <_topology->nbDomain(); idomain++) +// { +// MEDMEM::MESHING* mesh_builder=static_cast (_mesh[idomain]); +// if ( new_supports[idomain] ) +// mesh_builder->addGroup(*new_groups[idomain]); +// } +// //groups are copied by the addGroup method, +// //so they can be safely deleted here +// for (int i=0; i<_topology->nbDomain(); i++) +// { +// if ( new_supports[i] ) new_groups[i]->removeReference(); +// } + +// }// on groups +// }//on entities +// } + + +// void MESHCollection::castSupport(const MESHCollection& old_collection, vector& old_support, vector& new_support) +// { + +// if (old_collection._topology->nbDomain() != old_support.size()) +// { +// throw MED_EXCEPTION(STRING("Error : wrong call to MESHCollection::castSupport")); +// } +// vector > element_array (_topology->nbDomain()); + +// //parameters stored for passing description +// // from the old meshes to the new ones +// string name; +// string description; +// MED_EN::medEntityMesh entity; +// vector support_name(1); +// support_name[0]="support"; +// for (int inew=0; inew< _topology->nbDomain(); inew++) +// element_array[inew].clear(); + +// for (int idomain=0; idomain < old_collection._topology->nbDomain(); idomain++) +// { +// //retrieves the group igroup on domain idomain +// const MEDMEM::SUPPORT* support = old_support[idomain]; +// if (old_support[idomain]==0) continue; +// name = support->getName(); +// description=support->getDescription(); +// int nbelem = support->getNumberOfElements(MED_EN::MED_ALL_ELEMENTS); +// if (nbelem==0 && !_create_empty_groups) continue; + +// int* list_of_elems; +// if (support->isOnAllElements()) +// { +// list_of_elems = new int[nbelem]; +// for (int i=0; i (support->getNumber(MED_EN::MED_ALL_ELEMENTS)); + +// int* array=new int[nbelem]; +// int* ip=0; +// int* local=0; +// int* full_array=0; +// entity = support->getEntity(); +// int size; + +// switch (entity) +// { +// case MED_EN::MED_CELL : +// ip=new int[nbelem]; +// local= new int[nbelem]; +// size=nbelem; +// old_collection.getTopology()->convertCellToGlobal(idomain,list_of_elems,nbelem,array); +// _topology->convertGlobalCellList(array,nbelem,local,ip); +// for (int i=0; iconvertFaceToGlobal(idomain,list_of_elems,nbelem,array); +// _topology->convertGlobalFaceListWithTwins(array,nbelem,local,ip,full_array,size); +// for (int i=0; iconvertNodeToGlobal(idomain,list_of_elems,nbelem,array); +// _topology->convertGlobalNodeListWithTwins(array,nbelem,local,ip,full_array,size); +// for (int i=0; iisOnAllElements()) delete[] list_of_elems; +// } + +// //creating new groups from the previous list of elements +// for (int idomain=0; idomain <_topology->nbDomain(); idomain++) +// { +// if ( _mesh[idomain]->getNumberOfNodes() < 1 || +// (element_array[idomain].empty() && !_create_empty_groups)) +// { +// new_support[idomain]->removeReference(); +// new_support[idomain]=0; +// continue; +// } +// MEDMEM::SUPPORT* support= new_support[idomain]; +// support->setName(name); +// support->setMesh(_mesh[idomain]); +// support->setDescription(description); +// support->setEntity(entity); + +// element_array[idomain].sort(); +// element_array[idomain].unique(); + +// if ( element_array[idomain].empty() ) +// { +// support->setNumberOfGeometricType(0); +// } +// else +// { +// if (entity != MED_EN::MED_NODE) +// support->fillFromElementList(element_array[idomain]); +// else +// { +// support->fillFromNodeList(element_array[idomain]); +// } +// } +// } +// } + +// void MESHCollection::castField(const MESHCollection& old_collection, const string& fieldname, int itnumber, int ordernumber) +// { +// int type=old_collection.getDriver()->getFieldType(fieldname); +// char field_char[80]; +// strcpy(field_char,fieldname.c_str()); + +// if (type ==0) +// castFields(old_collection, field_char, itnumber, ordernumber); +// else +// castFields(old_collection, field_char, itnumber, ordernumber); +// } + +// void MESHCollection::castAllFields(const MESHCollection& initial_collection) +// { +// vector field_names; +// vector iternumber; +// vector ordernumber; +// vector types; +// initial_collection.getDriver()->readFileStruct(field_names,iternumber,ordernumber,types); + +// for (int i=0; i(initial_collection, field_char, iternumber[i], ordernumber[i]); +// else +// castFields(initial_collection, field_char, iternumber[i], ordernumber[i]); +// } +// } + +// //! \param initial_collection mesh collection before the redistribution stage +// //! \param idomain domain number on the new collection +// //! \param entity dimension of the treated nodes (0 : cell, 1 : face in 3d, edge in 2d, 2: edge in 3d) +// void MESHCollection::createNodalConnectivity(const MESHCollection& initial_collection,int idomain, int entity) +// { + + +// int dimension=0; +// int nb_elems=0; +// ParaMEDMEM::MEDCouplingUMesh* mesh_builder = m_mesh[idomain]; + +// //creating arrays for storing global numbers and cell types +// switch (entity) +// { +// case MED_EN::MED_CELL: +// dimension=initial_collection.getMeshDimension(); +// nb_elems=m_topology->getCellNumber(idomain); +// break; +// case MED_EN::MED_EDGE: +// case MED_EN::MED_FACE: +// dimension=initial_collection.getMeshDimension()-1; +// nb_elems=m_topology->getFaceNumber(idomain); +// break; +// default: +// nb_elems=0; +// break; +// } + +// if (nb_elems == 0) return; +// SCRUTE_MED(nb_elems); + + +// int *list= new int[nb_elems]; +// MED_EN::medGeometryElement *cell_type_list= new MED_EN::medGeometryElement[nb_elems]; + + +// // cout << "Beginning of retrieval "<getCellList(idomain,list); +// break; +// case MED_EN::MED_EDGE: +// case MED_EN::MED_FACE: +// m_topology->getFaceList(idomain,list); +// break; +// default: + +// break; +// } + +// //retrieving cell_types +// initial_collection.getTypeList(list,nb_elems,entity,cell_type_list); +// // cout <<"end of type retrieval"<::iterator iter= type_numbers.find(cell_type_list[icell]); +// if (iter!=type_numbers.end()) +// (iter->second)++; +// else +// type_numbers[cell_type_list[icell]]=1; + +// } +// //cout << "Nombre de tetras"< index; + +// map type_cell_list; + +// MED_EN::MESH_ENTITIES::const_iterator currentEntity; +// std::map::const_iterator iter; +// //currentEntity = MED_EN::meshEntities.find(entity); +// for (iter = type_numbers.begin();iter != type_numbers.end(); iter++) +// { +// MED_EN::medGeometryElement type = iter->first; +// if (!isDimensionOK(type,dimension)) continue; +// //if (iter->second==0) continue; +// index[type]=0; +// type_cell_list[type]=new int[type_numbers[type]]; +// // cout << "type :"<first; +// if (!isDimensionOK(type,dimension)) continue; + +// type_array[itype]=type; + +// present_type_numbers[itype]=type_numbers[type]; + +// MESSAGE_MED("Nombre d'elements de type "<get +// map type_connectivity; +// vector polygon_conn; +// vector polygon_conn_index; +// vector polyhedron_conn; +// vector polyhedron_conn_index; +// vector polyhedron_face_index; + +// //Treating nodes + +// DataArrayInt* conn_index = (initial_collection.getMesh())[idomain].getConnIndex(); +// DataArrayInt* index=(initial_collection.getMesh())[idomain].getIndex(); + + +// _topology->convertGlobalCellList(cell_list,nb_cells,local,ip); +// for (iter = type_numbers.begin();iter != type_numbers.end(); iter++) +// { +// MED_EN::medGeometryElement type = iter->first; + + +// if (!isDimensionOK(type,dimension)) continue; +// //if (type_numbers[type]==0) continue; +// if (type != MED_EN::MED_POLYGON && type != MED_EN::MED_POLYHEDRA) +// { +// int nbnode_per_type = (int)type%100; +// type_connectivity[type]=new int[type_numbers[type]*nbnode_per_type]; +// initial_collection.getNodeConnectivity(type_cell_list[type],type_numbers[type],entity,type,type_connectivity[type]); +// } +// else if (type == MED_EN::MED_POLYGON && dimension==2) +// { +// initial_collection.getPolygonNodeConnectivity(type_cell_list[type],type_numbers[type],entity,polygon_conn,polygon_conn_index); +// //type_connectivity[type]=&polygon_conn[0]; +// } +// else if (type == MED_EN::MED_POLYHEDRA && dimension==3) +// { +// initial_collection.getPolyhedraNodeConnectivity(type_cell_list[type],type_numbers[type],entity,polyhedron_conn, +// polyhedron_conn_index, polyhedron_face_index); +// //type_connectivity[type]=&polygon_conn[0]; +// } +// delete[] type_cell_list[type]; +// } + +// //creating node mapping +// //!TODO : compute the total number of nodes +// if (entity==MED_EN::MED_CELL) +// { +// m_topology->createNodeMapping(type_connectivity,type_numbers,polygon_conn,polygon_conn_index, +// polyhedron_conn,polyhedron_conn_index,polyhedron_face_index,idomain); +// } + +// //converting node global numberings to local numberings +// //for (iter = (*currentEntity).second.begin();iter != (*currentEntity).second.end(); iter++) +// for (iter = type_numbers.begin();iter != type_numbers.end(); iter++) +// { +// MED_EN::medGeometryElement type = iter->first; + +// if (!isDimensionOK(type, dimension)) continue; +// if (type_numbers[type]==0) continue; +// if (type != MED_EN::MED_POLYGON && type != MED_EN::MED_POLYHEDRA) +// { +// int nbnode_per_type = (int)type%100; +// m_topology->convertToLocal2ndVersion(type_connectivity[type],type_numbers[type]*nbnode_per_type,idomain); +// } +// else if (type == MED_EN::MED_POLYGON && dimension==2) +// { +// int nbpoly = type_numbers[type]; +// m_topology->convertToLocal2ndVersion(&polygon_conn[0], polygon_conn_index[nbpoly]-1, idomain); +// } +// else if (type == MED_EN::MED_POLYHEDRA && dimension==3) +// { +// int nbpoly = type_numbers[type]; +// m_topology->convertToLocal2ndVersion(&polyhedron_conn[0], polyhedron_face_index[polyhedron_conn_index[nbpoly]-1]-1, idomain); +// } + +// } + + +// //writing coordinates +// if (entity==MED_EN::MED_CELL) +// { +// //setting coordinates from initial_collection coordinates +// int nbnode=m_topology->getNodeNumber(idomain); +// MESSAGE_MED("Number of nodes on domain "<< idomain <<" : "<getNodeList(idomain,node_list); +// initial_collection.getCoordinates(node_list,nbnode,coordinates); +// delete[] node_list; + +// // redundant specification of number of nodes is required!! MED imperfection, sorry... + +// mesh_builder->setNumberOfNodes(nbnode); +// //TODO : change MEDMEM so that it accepts a direct setting of coordinates +// // (in the present version, it is deep-copied) +// mesh_builder->setCoordinates(initial_collection.getSpaceDimension(), +// m_topology->getNodeNumber(idomain), coordinates, initial_collection.getSystem(), +// MED_EN::MED_FULL_INTERLACE); +// delete [] coordinates; +// } + +// int nb_plain_types=0; +// for (iter = type_numbers.begin();iter != type_numbers.end(); iter++) +// { +// MED_EN::medGeometryElement type = iter->first; + +// if (!isDimensionOK(type, dimension)) continue; +// if (type_numbers[type]==0) continue; +// if (type != MED_EN::MED_POLYGON && type != MED_EN::MED_POLYHEDRA) +// nb_plain_types++; +// } +// mesh_builder->setNumberOfTypes(nb_plain_types,entity); +// mesh_builder->setTypes(type_array,entity); +// mesh_builder->setNumberOfElements(present_type_numbers,entity); +// if (entity==MED_EN::MED_CELL) +// mesh_builder->setMeshDimension(dimension); + +// delete[]present_type_numbers; +// delete[]type_array; +// //setting node connectivities +// for (iter = type_numbers.begin();iter != type_numbers.end(); iter++) +// { +// MED_EN::medGeometryElement type = iter->first; + +// if (!isDimensionOK(type,dimension)) continue; +// if (type_numbers[type]==0) continue; + +// if (type != MED_EN::MED_POLYHEDRA && type != MED_EN::MED_POLYGON) +// { +// mesh_builder->setConnectivity(type_connectivity[type],entity,type); +// delete[] type_connectivity[type]; +// } +// else if (type == MED_EN::MED_POLYGON && dimension ==2) +// { +// mesh_builder->setPolygonsConnectivity(&polygon_conn_index[0], +// &polygon_conn[0], +// type_numbers[type], +// entity); +// } +// else if (type == MED_EN::MED_POLYHEDRA && dimension ==3) +// { +// mesh_builder->setPolyhedraConnectivity(&polyhedron_conn_index[0], +// &polyhedron_face_index[0], +// &polyhedron_conn[0], +// type_numbers[type], +// entity); + +// } +// } +// MESSAGE_MED("end of createNodalConnectivity"); +// } + + +// /*! retrieves the faces that are present in a mesh and stores them in a +// * dynamic structure made of a map of MEDPARTITIONER_FaceModel +// * +// * \param idomain domain id on which the faces are collected +// * \param face_map container storing the faces +// */ +// void MESHCollection::getFaces(int idomain, +// map >& face_map) +// { +// MED_EN::medEntityMesh constituent_entity = getSubEntity(); +// const medGeometryElement* types; +// try +// { +// types = _mesh[idomain]->getTypes(constituent_entity); +// } +// catch(MEDEXCEPTION&){ return;} + +// int nbtypes = _mesh[idomain]->getNumberOfTypes(constituent_entity); +// const int* global_numbering= _mesh[idomain]->getGlobalNumberingIndex(constituent_entity); +// int* conn = const_cast (_mesh[idomain]->getConnectivity(MED_EN::MED_FULL_INTERLACE,MED_EN::MED_NODAL,constituent_entity, MED_EN::MED_ALL_ELEMENTS)); +// for (int itype=0; itypesetType(type); +// int nbnodes = type%100; +// face_model->setNbNodes(nbnodes); +// face_model->setGlobal(_topology->convertFaceToGlobal(idomain,iface)); +// for (int i=0; igetElementType(MED_EN::MED_CELL,ilocal1); +// MEDMEM::CELLMODEL celltype1 (type1); + +// const int* conn_index1 = _mesh[ip1]->getConnectivityIndex(MED_EN::MED_NODAL,MED_EN::MED_CELL); +// const int* conn1 = _mesh[ip1]->getConnectivity(MED_EN::MED_FULL_INTERLACE,MED_EN::MED_NODAL,MED_EN::MED_CELL,MED_EN::MED_ALL_ELEMENTS); + +// // MED_EN::medGeometryElement type2 = _mesh[ip2]->getElementType(MED_EN::MED_CELL,ilocal2); +// //MEDMEM::CELLTYPE celltype2 (type2); +// const int* conn_index2 = _mesh[ip2]->getConnectivityIndex(MED_EN::MED_NODAL,MED_EN::MED_CELL); +// const int* conn2 = _mesh[ip2]->getConnectivity(MED_EN::MED_FULL_INTERLACE,MED_EN::MED_NODAL,MED_EN::MED_CELL,MED_EN::MED_ALL_ELEMENTS); + +// vector nodes1, nodes1_local; +// vector nodes2; +// for (int i= conn_index1[ilocal1-1]; iconvertNodeToGlobal(ip1,*(conn1+i-1))); +// nodes1_local.push_back( conn1[i-1] ); +// } +// for (int i= conn_index2[ilocal2-1]; iconvertNodeToGlobal(ip2,*(conn2+i-1))); + +// return MEDPARTITIONER_FaceModel::getCommonFace( &nodes1[0], &nodes1_local[0], celltype1, +// &nodes2[0], nodes2.size(), face_index); +// } + +// //================================================================================ +// /*! +// * \brief Makes a face common for two given cells +// * \param nodes1 - globl nodes of the first cell +// * \param nodes1_local - local nodes of the first cell +// * \param celltype1 - cell model of the first cell +// * \param nodes2 - globl nodes of the second cell +// * \param nb_nodes2 - nb of nodes of the second cell +// * \param global_id - id of the new common face +// */ +// //================================================================================ + +// MEDPARTITIONER_FaceModel* +// MEDPARTITIONER_FaceModel::getCommonFace(const int* nodes1, +// const int* nodes1_local, +// const MEDMEM::CELLMODEL& celltype1, +// const int* nodes2, +// int nb_nodes2, +// int global_id) +// { +// int nbfaces= celltype1.getNumberOfConstituents(1); +// int ** faces = celltype1.getConstituents(1); +// MED_EN::medGeometryElement* types = celltype1.getConstituentsType(1); +// int iface=0; +// int dimension=celltype1.getDimension(); + +// while (iface=dimension) break; +// iface++; +// } + +// if (iface==nbfaces) +// throw MEDEXCEPTION("MEDPARTITIONER::getCommonFace - No common face found !"); + +// MEDPARTITIONER_FaceModel* face_model = new MEDPARTITIONER_FaceModel; +// face_model->setType(types[iface]); +// int nbnodes = types[iface]%100; +// face_model->setNbNodes(nbnodes); +// face_model->setGlobal(global_id); +// for (int i=0; i +#include + +namespace ParaMEDMEM +{ + class MEDCouplingUMesh; +} + +namespace MEDPARTITIONER +{ + + class Topology; + class MESHCollectionDriver; + class ParaDomainSelector; + class MEDSKYLINEARRAY; + class CONNECTZONE; + + typedef enum{MedAscii, MedXML, Undefined} DriverType; + + + class MEDPARTITIONER_EXPORT MESHCollection + { + + public: + + //Default constructor + MESHCollection(); + + //Constructing from an existing mesh and a new topology + MESHCollection(MESHCollection&, Topology*, bool family_splitting=false, bool create_empty_groups=false); + + //Constructing the mesh collection from a file + MESHCollection(const std::string& filename); + + //Constructing the mesh collection from a file + MESHCollection(const std::string& filename, ParaDomainSelector& domainSelector); + + //Constructing the mesh collection from a file + MESHCollection(const std::string& filename, const std::string& meshname); + + ~MESHCollection(); + + bool isParallelMode() const { return _domain_selector; } + + //writing to a distributed file + void write(const std::string& filename); + + //getting the driver + MESHCollectionDriver* retrieveDriver(); + MESHCollectionDriver* getDriver() const; + void setDriverType(MEDPARTITIONER::DriverType type) {_driver_type=type;} + + //creation of the cell graph + void buildCellGraph(MEDPARTITIONER::MEDSKYLINEARRAY* & array,int *& edgeweights ); + + //creation and partition of the associated graph + Topology* createPartition(int nbdomain, Graph::splitter_type type = Graph::METIS, + const std::string& ="", int* edgeweights=0, int* verticesweights=0); + + //creation of a user specified partition + Topology* createPartition(const int* partition); + + // //retrieving list of types +// void getTypeList(int* cell_list,int nb_cells,MED_EN::medEntityMesh entity, +// MED_EN::medGeometryElement* type_list) const ; + +// //getting list of coordinates +// void getCoordinates(int* node_list,int nb_nodes, double* coordinates) const ; + +// //getting connectivities +// void getNodeConnectivity( const int* cell_list,int nb_cells,MED_EN::medEntityMesh,MED_EN::medGeometryElement type, int* type_connectivity) const ; +// void getPolygonNodeConnectivity(const int* cell_list,int nb_cells,MED_EN::medEntityMesh entity, +// vector& type_connectivity, vector& connectivity_index) const; +// void getPolyhedraNodeConnectivity(const int* cell_list,int nb_cells,MED_EN::medEntityMesh entity, +// vector& type_connectivity, vector& connectivity_index, vector& face_connectivity_index) const; + +// void getFaceConnectivity( const int* cell_list,int nb_cells,MED_EN::medEntityMesh,MED_EN::medGeometryElement type, int* type_connectivity) const ; + +// //void getFaceConnectivity( const int* cell_list,int nb_cells,MED_EN::medGeometryElement type, int* type_connectivity) const ; + +// //getting mesh dimension + int getMeshDimension() const ; + +// //getting space dimension +// int getSpaceDimension() const ; + +// //getting system of coordinates +// std::string getSystem() const; + +// //getting name of the mesh +// std::string getMeshName() const; + +// //return constituent entity +// MED_EN::medEntityMesh getSubEntity() const; + +// //getting a reference to mesh vector + std::vector& getMesh() ; + std::vector& getFaceMesh() ; + std::vector >&getGroupMeshes(); + + ParaMEDMEM::MEDCouplingUMesh* getMesh(int idomain); + ParaMEDMEM::MEDCouplingUMesh* getFaceMesh(int idomain); + std::vector& getGroupMeshes(int idomain); + +// //getting a reference to a domain mesh +// MEDMEM::MESH* getMesh(int) const; + +// //getting a reference to connect zones vector + std::vector& getCZ(); + + //getting a pointer to topology + Topology* getTopology() const ; + + + //settig a new topology + void setTopology(Topology* topology); + + //getting/setting the name of the global mesh (as opposed + //to the name of a subdomain \a nn, which is name_nn) + std::string getName() const {return _name;} + void setName(const std::string& name) {_name=name;} + + //getting/setting the description of the global mesh + std::string getDescription() const {return _description;} + void setDescription(const std::string& name) { _description=name;} + + //creates the node mapping between an old collection and the present one + void createNodeMapping( MESHCollection& initialCollection, std::multimap,std::pair >& nodeMapping); + + //creates faces on the new collection + void castFaceMeshes( MESHCollection& initialCollection,const std::multimap,std::pair >& nodeMapping); + + //creates faces on the new collection + void castMeshes(std::vector& meshCastFrom,std::vector& meshCastTo, MESHCollection& initialCollection,const std::multimap,std::pair >& nodeMapping); + + //!transfers families from an old MESHCollection to new mesh + // void castFamilies(const MESHCollection& old_collection); + + //void castSupport(const MESHCollection& old_collection, std::vector& old_support, std::vector& new_support); + +// //!casts all the fields to the new mesh collection +// void castAllFields(const MESHCollection& old_collection); + +// //!casts one specific field to the new mesh collection +// void castField(const MESHCollection& old_collection, const string& fieldname, int itnumber, int ordernumber); + +// //choosing a group to be indivisible +// void setIndivisibleGroup(const string& a); + +// //!constructing connect zones +// // void buildConnectZones(int idomain); +// void buildConnectZones(); +// void buildConnectZonesBetweenProcs(std::vector > > &, map< pair, MEDPARTITIONER::MEDSKYLINEARRAY*> & local_cell_cell_correspondency); + +// void addJointGroup(const std::vector&, int idomain, int idistant); + +// static bool isDimensionOK(MED_EN::medGeometryElement type, int dim) +// { +// return ((type/100 == dim) || (dim==2 && type == MED_EN::MED_POLYGON) || (dim==3 && type == MED_EN::MED_POLYHEDRA)); +// } +// void setSubdomainBoundaryCreates(bool flag) { _subdomain_boundary_creates=flag;} +// bool getSubdomainBoundaryCreates(){return _subdomain_boundary_creates;} + +// void setFamilySplitting(bool flag){_family_splitting=flag;} +// bool getFamilySplitting(){return _family_splitting;} + +// void setCreateEmptyGroups(bool flag){_create_empty_groups=flag;} +// bool getCreateEmptyGroups(){return _create_empty_groups;} + + + private: + + // //!creates connectivities for a domain and an entity (face or cell) +// void createNodalConnectivity(const MESHCollection & initial_collection, int idomain, MED_EN::medEntityMesh entity); + +// //!creates the tags for indivisible groups +// void treatIndivisibleRegions(int* tag); + +// //!projects a field from an old collection to the present one +// //!field is identified by (name, dt, it) +// template +// void castFields(const MESHCollection& old_collection, const string& fieldname, int itnumber, int ordernumber); + +// void getFaces(int idomain, std::map >&); + +// MEDPARTITIONER_FaceModel* getCommonFace(int ip1,int ilocal1,int ip2,int ilocal2,int face_index); + +// template +// void fillGlobalConnectivity(TID2CONN & node2cell, TID2CONN & cell2node ); + + //!link to mesh_collection topology + Topology* _topology; + + //!control over topology + bool _owns_topology; + + //!link to graph + // Graph* _cell_graph; + + //! Driver for read/write operations + MESHCollectionDriver* _driver; + + //! Parallelizer - mark of parallel execution mode + ParaDomainSelector* _domain_selector; + + //!links to meshes + std::vector _mesh; + std::vector _faceMesh; + std::vector< std::vector > _groupMesh; + + //!index of a non empty mesh within _mesh (in parallel mode all of meshes can be empty) + int _i_non_empty_mesh; + + //!links to connectzones + std::vector _connect_zones; + + //!list of groups that are not to be splitted + std::vector _indivisible_regions; + + //!name of global mesh + std::string _name; + + //!description of global mesh + std::string _description; + + //! specifies the driver associated to the collection + DriverType _driver_type; + + /*! flag specifying that the splitter should create boundary constituent entity + so that they are written in joints*/ + bool _subdomain_boundary_creates; + + /*! flag specifying that families must be preserved by the + splitting*/ + bool _family_splitting; + + /*! flag specifying that groups must be created on all domains, + even if they are empty*/ + bool _create_empty_groups; + }; + +}//of namespace + +#endif /*MESHCOLLECTION_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.cxx b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.cxx new file mode 100644 index 000000000..3a6aa35d1 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.cxx @@ -0,0 +1,697 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#include +#include +#include +#include + +#include +#include + +#include +#include +#include +#include + +#include +#include "MEDCouplingUMesh.hxx" +#include "MEDLoader.hxx" +#include "MEDMEM_Exception.hxx" +extern "C" { +#include "med.h" +} +//MEDPARTITIONER includes +#include "MEDPARTITIONER_ConnectZone.hxx" +#include "MEDPARTITIONER_Topology.hxx" +#include "MEDPARTITIONER_ParallelTopology.hxx" +#include "MEDPARTITIONER_SequentialTopology.hxx" +#include "MEDPARTITIONER_MESHCollectionDriver.hxx" +#include "MEDPARTITIONER_MESHCollection.hxx" +#include "MEDPARTITIONER_ParaDomainSelector.hxx" + +using namespace MEDPARTITIONER; + +//template inclusion +//#include "MEDPARTITIONER_MESHCollectionDriver.H" + +med_geometrie_element typmai[MED_NBR_GEOMETRIE_MAILLE+2] = { MED_POINT1, + MED_SEG2, + MED_SEG3, + MED_TRIA3, + MED_TRIA6, + MED_QUAD4, + MED_QUAD8, + MED_TETRA4, + MED_TETRA10, + MED_HEXA8, + MED_HEXA20, + MED_PENTA6, + MED_PENTA15, + MED_PYRA5, + MED_PYRA13, + MED_POLYGONE, + MED_POLYEDRE }; + +MESHCollectionDriver::MESHCollectionDriver(MESHCollection* collection):_collection(collection) +{ +} + + +/*!reads a unique MED File v>=2.1 + * and mounts the corresponding mesh in memory + *\param filename binary file + *\param meshname mesh name in the MED file + * */ +int MESHCollectionDriver::readSeq(char* filename, char* meshname) +{ + //const char* LOC = "MEDPARTITIONER::MESHCollectionDriver::readSeq()"; + // BEGIN_OF_MED(LOC); + + _filename.resize(1); + _filename[0]=string(filename); + //puts the only mesh in the mesh vector + (_collection->getMesh()).push_back(MEDLoader::ReadUMeshFromFile(filename, meshname, 0)); + (_collection->getFaceMesh()).push_back(MEDLoader::ReadUMeshFromFile(filename,meshname,-1)); + + //reading groups + vector groupNames = MEDLoader::GetMeshGroupsNames(filename,meshname); + (_collection->getGroupMeshes()).resize(groupNames.size()); + + for (int i=0; i< groupNames.size();i++) + { + vector myGroup; + myGroup.push_back(groupNames[i]); + (_collection->getGroupMeshes())[i].push_back(MEDLoader::ReadUMeshFromGroups(filename,meshname,-1,myGroup)); + } + + + _collection->setName(meshname); + (_collection->getCZ()).clear(); + vector cellglobal,nodeglobal,faceglobal; + cellglobal.resize(1); + nodeglobal.resize(1); + faceglobal.resize(1); + cellglobal[0]=0; + nodeglobal[0]=0; + faceglobal[0]=0; + //creation of topology from mesh + //connectzone argument is 0 + ParallelTopology* aPT = new ParallelTopology + ((_collection->getMesh()), (_collection->getCZ()), cellglobal, nodeglobal, faceglobal); + _collection->setTopology(aPT); + // END_OF_MED(LOC); + return 0; +} + +// /*! +// * Reads the file structure to determine the list +// * of all the available fields +// * +// * \param field_names, vector containing the field names +// * \param iternumber, vector containing the iteration numbers +// * \param ordernumber, vector containing the order numbers +// * \param types, vector containing 0 for int fields and 1 for double fields +// * +// */ + +// void MESHCollectionDriver::readFileStruct(vector & field_names,vector& iternumber,vector & ordernumber, vector & types) +// { +// const char* LOC = "MEDPARTITIONER::MESHCollectionDriver::readFileStruct()"; +// BEGIN_OF_MED(LOC); + +// const MEDMEM::MED med_struct (MEDMEM::MED_DRIVER,_filename[0]); +// int nb_fields = med_struct.getNumberOfFields(); + +// MESSAGE_MED("found "< names = med_struct.getFieldNames(); +// for (int ifield = 0; ifield < nb_fields; ifield++) +// { +// deque dtit=med_struct.getFieldIteration(names[ifield]); + +// for (deque::const_iterator iter =dtit.begin(); iter!=dtit.end(); iter++) +// { +// field_names.push_back(names[ifield]); +// iternumber.push_back(iter->dt); +// ordernumber.push_back(iter->it); + +// // testing whether the field is of double or int type +// MEDMEM::FIELD_* field = med_struct.getField(names[ifield],iter->dt,iter->it); +// if (dynamic_cast*>(field)) +// types.push_back(1); +// else +// types.push_back(0); + +// } +// } +// END_OF_MED(LOC); +// } + +// //!retrieves the type of a field for a given fieldname +// int MESHCollectionDriver::getFieldType(const string& fieldname) +// { +// const char* LOC = "MEDPARTITIONER::MESHCollectionDriver::getFieldType()"; +// BEGIN_OF_MED(LOC); +// const MEDMEM::MED med_struct (MEDMEM::MED_DRIVER,_filename[0]); + +// deque dtit=med_struct.getFieldIteration(fieldname); +// deque::const_iterator iter =dtit.begin(); + +// // testing whether the field is of double or int type +// MEDMEM::FIELD_* field = med_struct.getField(fieldname,iter->dt,iter->it); + +// END_OF_MED(LOC); + +// if (dynamic_cast*>(field)) +// return 1; +// else +// return 0; + +// } +// namespace +// { +// //================================================================================ +// /*! +// * \brief Structure used in the method below +// */ +// struct TJointData +// { +// char _name[MED_TAILLE_NOM+1]; +// int _nb_couples; +// med_int _distant_domain; +// med_geometrie_element _geo_local, _geo_dist; +// }; +// } + +// //================================================================================ +// /*! +// * \brief Read CELL-CELL correspondences of joints with domains on other procs +// * \param idomain - domain index to return correspondence for +// * \param loc_domains - domians on this pocessor +// * \param domain_selector - info on cell distribution among procs +// * \param loc2glob_corr - out, correspondence pairs where distant ids are global +// */ +// //================================================================================ + +// void MESHCollectionDriver::readLoc2GlobCellConnect(int idomain, +// const set& loc_domains, +// ParaDomainSelector* domain_selector, +// vector & loc2glob_corr) +// { +// using namespace med_2_3; +// med_err err; + +// // find joints of domains touching idomain and loaded on other processors + +// TJointData joint; +// list< TJointData > joints; +// int total_nb_couples = 0; + +// MEDMEM::MESH* loc_mesh = _collection->getMesh()[idomain]; +// char* meshname = (char*) _meshname[idomain].c_str(); +// char* filename = (char*) _filename[idomain].c_str(); +// //cout << "#" << domain_selector->rank() << ": mesh - " << meshname << endl; + +// med_idt fid = MEDouvrir( filename, MED_LECTURE); +// int njoint = MEDnJoint(fid, meshname); +// for (int ijoint=0; ijointgetTypesWithPoly(MED_EN::MED_CELL); +// int nbtypes = loc_mesh->getNumberOfTypesWithPoly(MED_EN::MED_CELL); +// const list& all_types = MED_EN::meshEntities[ MED_EN::MED_CELL ]; +// for (int itype=0; itype::const_iterator dist_type = all_types.begin(); +// for ( ; dist_type != all_types.end(); ++dist_type ) +// { +// if ( !_collection->isDimensionOK( *dist_type, loc_mesh->getMeshDimension() )) +// continue; +// joint._geo_dist = (med_geometrie_element) *dist_type; +// joint._nb_couples = MEDjointnCorres(fid, meshname, joint._name, +// MED_MAILLE, joint._geo_local, +// MED_MAILLE, joint._geo_dist ); +// if ( joint._nb_couples > 0 ) +// { +// joints.push_back( joint ); +// total_nb_couples += joint._nb_couples; +// } +// } +// } +// delete [] types; +// } + +// // read cell pairs of found joints and replace distant local ids with global ones + +// loc2glob_corr.resize( 2 * total_nb_couples ); +// if ( total_nb_couples > 0 ) +// { +// int* cell_corresp = & loc2glob_corr[0]; + +// list< TJointData >::iterator jnt = joints.begin(); +// for ( ; jnt != joints.end(); ++jnt ) +// { +// // read cell couples +// if ( MEDjointLire( fid, meshname, jnt->_name, +// cell_corresp, 2 * jnt->_nb_couples, +// MED_MAILLE, jnt->_geo_local, +// MED_MAILLE, jnt->_geo_dist ) < 0 ) continue; + +// // distant local ids -> global ids +// if ( int shift_to_global = domain_selector->getDomainShift( jnt->_distant_domain )) +// for ( int i = 0 ; i < jnt->_nb_couples; ++i ) +// cell_corresp[ 2*i + 1 ] += shift_to_global; + +// cell_corresp += 2 * jnt->_nb_couples; +// } +// } +// MEDfermer( fid ); +// } + +//================================================================================ +/*! + * \brief Return mesh dimension from distributed med file had being read + */ +//================================================================================ + +// int MESHCollectionDriver::readMeshDimension() const +// { +// const char* LOC = "MESHCollectionDriver::readMeshDimension(): "; +// if ( _filename.empty() || _meshname.empty() ) +// throw MEDMEM::MED_EXCEPTION( "file name or mesh name not available"); + +// MEDMEM::MED med(MEDMEM::MED_DRIVER, _filename[0]); +// if ( MEDMEM::MESH* mesh = med.getMesh( _meshname[0] )) +// return mesh->getMeshDimension(); + +// throw MEDMEM::MED_EXCEPTION( MEDMEM::STRING(LOC) << "mesh name is invalid"); +// } + +void MESHCollectionDriver::readSubdomain(vector& cellglobal, + vector& faceglobal, + vector& nodeglobal, int idomain) +{ +// const char* LOC = "MEDPARTITIONER::MESHCollectionDriver::readSubdomain()"; +// BEGIN_OF_MED(LOC); + char file[256]; + char meshname[MED_TAILLE_NOM+1]; + + strcpy(meshname,_meshname[idomain].c_str()); + strcpy(file,_filename[idomain].c_str()); + cout << "Reading "<<_meshname[idomain]<<" in "<<_filename[idomain]<getMesh())[idomain]=new MEDMEM::MESH(MEDMEM::MED_DRIVER,file, meshname); + (_collection->getMesh())[idomain]=MEDLoader::ReadUMeshFromFile(file,meshname,0); + (_collection->getFaceMesh())[idomain]=MEDLoader::ReadUMeshFromFile(file,meshname,-1); + + //reading groups + vector groupNames = MEDLoader::GetMeshGroupsNames(file,meshname); + if (idomain==0) + (_collection->getGroupMeshes()).resize(groupNames.size()); + + for (int i=0; i< groupNames.size();i++) + { + vector myGroup; + myGroup.push_back(groupNames[i]); + (_collection->getGroupMeshes())[i].push_back(MEDLoader::ReadUMeshFromGroups(file,meshname,-1,myGroup)); + } + + cout <<"End of Read"<setName(string(name)); + cz->setDescription(joint_description); + cz->setLocalDomainNumber(idomain); + cz->setDistantDomainNumber(distant); + cz->setLocalMesh((_collection->getMesh())[idomain]); + cz->setDistantMesh((_collection->getMesh())[distant]); + cz->setNodeCorresp(node_corresp,ncouples); + (_collection->getCZ()).push_back(cz); + + }//loop on correspom_topology->nbDomain())ndances + }//loop on joints + + // Reading global numbering + // à faire dans MEDLoader + + // MEDCouplingFieldDouble globalcell=MEDLoader::GetGlobalNumbering(file,meshname,0); + + + MEDfermer(fid); + + // END_OF_MED(LOC); +} + +void MESHCollectionDriver::writeSubdomain(int idomain, int nbdomains, char* distfilename, + ParaDomainSelector* domain_selector) +{ + //build connect zones + // if (nbdomains>1) + // _collection->buildConnectZones(idomain); + + // MESSAGE_MED(" Number of connect zones "<<(_collection->getCZ()).size()); + + //writing connect zones in joints + + med_idt fid = MEDouvrir(distfilename,MED_LECTURE_ECRITURE); + + int index_joint=0; + + + for (int icz=0; icz<(_collection->getCZ()).size(); icz++) + { + if ((_collection->getCZ())[icz]->getLocalDomainNumber()==idomain) + { + med_err error; + int idistant=(_collection->getCZ())[icz]->getDistantDomainNumber(); + char joint_name[MED_TAILLE_NOM+1]; + sprintf(joint_name,"joint_%i",idistant+1); + char desc[MED_TAILLE_DESC+1]; + sprintf(desc,"connect_zone_%d",icz+1); + + char distant_name[MED_TAILLE_NOM+1]; + //sprintf(distant_name,"domain_%i",(_collection->getCZ())[icz]->getDistantDomainNumber()); + + // sprintf(distant_name,(_collection->getMesh())[idistant]->getName().c_str()); + sprintf(distant_name,"domain_%i",idistant); + char mesh_name[MED_TAILLE_NOM+1]; + + strcpy (mesh_name, (_collection->getMesh())[idomain]->getName()); + // SCRUTE_MED((_collection->getMesh())[idomain]->getName()); + error = MEDjointCr(fid,mesh_name, joint_name, desc, + idistant, distant_name); + if (error==-1) cout << "erreur creation de joint "<getCZ())[icz]->getNodeNumber(); + int* node_corresp=const_cast((_collection->getCZ())[icz]->getNodeCorrespValue()); + + /* Nodes are reordered so that the ordering on the local and the distant domain + correspond. The chosen order is the natural ordering on the domain + with lowest proc id*/ + // if (_collection->getSubdomainBoundaryCreates()) +// if (idomaingetSubdomainBoundaryCreates()) +// { +// MED_EN::medEntityMesh constituent_entity = +// (_collection->getMeshDimension()==3)?MED_EN::MED_FACE:MED_EN::MED_EDGE; +// writeElementJoint(constituent_entity, icz, idomain, idistant, mesh_name,joint_name,fid); +// } + index_joint++; + } + } + + char meshchar[MED_TAILLE_NOM+1]; + strcpy(meshchar,(_collection->getMesh())[idomain]->getName()); + + // Writing cell global numbering +// // +// { +// int ncell=_collection->getTopology()->getCellNumber(idomain); +// int* array=new int[ncell]; +// _collection->getTopology()->getCellList(idomain,array); +// int offset=0; + + // MED_EN::MESH_ENTITIES::const_iterator currentEntity; +// list::const_iterator iter; +// currentEntity = MED_EN::meshEntities.find(MED_EN::MED_CELL); + +// int nbtypes = (_collection->getMesh())[idomain]->getNumberOfTypesWithPoly(MED_EN::MED_CELL); +// MED_EN::medGeometryElement* types =(_collection->getMesh())[idomain]->getTypesWithPoly(MED_EN::MED_CELL); +// for (int itype=0; itypeisDimensionOK(type,_collection->getMeshDimension())) continue; +// int ntype = (_collection->getMesh())[idomain]->getNumberOfElementsWithPoly(MED_EN::MED_CELL,type); +// if (ntype==0) continue; +// MEDglobalNumEcr(fid,meshchar, array+offset, ntype, +// MED_MAILLE, (med_geometrie_element)type); +// offset+=ntype; + +// } +// delete[] types; + // delete[] array; + // } + +// MED_EN::medEntityMesh constituent_entity; +// if (_collection->getMeshDimension()==3) +// constituent_entity=MED_EN::MED_FACE; +// else if (_collection->getMeshDimension()==2) +// constituent_entity=MED_EN::MED_EDGE; +// else throw MEDEXCEPTION("Wrong dimension"); + + +// //writing face global numbering +// { +// int * array; +// int offset=0; +// int nface= _collection->getTopology()->getFaceNumber(idomain); +// if (nface >0) +// array=new int[nface]; +// _collection->getTopology()->getFaceList(idomain,array); + +// int nbfacetypes = (_collection->getMesh())[idomain]->getNumberOfTypesWithPoly(constituent_entity); +// MED_EN::medGeometryElement* facetypes; +// if (nbfacetypes>0) +// facetypes =(_collection->getMesh())[idomain]->getTypesWithPoly(constituent_entity); + +// for (int itype=0; itypeisDimensionOK(type,_collection->getMeshDimension()-1)) continue; + +// int ntype = (_collection->getMesh())[idomain]->getNumberOfElementsWithPoly(constituent_entity,type); +// if (ntype==0) continue; +// MEDglobalNumEcr(fid,meshchar, array+offset, ntype, +// MED_MAILLE, (med_geometrie_element)type); + +// offset+=ntype; +// } +// if (nface>0) delete[] array; +// if (nbfacetypes>0) delete[] facetypes; +// } + + +// //writing node global numbering +// { +// int nnode=_collection->getTopology()->getNodeNumber(idomain); +// int* array = new int[nnode]; +// _collection->getTopology()->getNodeList(idomain,array); +// MEDglobalNumEcr(fid,meshchar, array, nnode, +// MED_NOEUD, MED_POINT1); +// delete[] array; +// } + + MEDfermer(fid); + std::cout<<"End of writing"<getCZ())[icz]->getEntityCorrespNumber(entity,entity); +// const int* index = (_collection->getCZ())[icz]->getEntityCorrespIndex(entity,entity); +// const int* value = (_collection->getCZ())[icz]->getEntityCorrespValue(entity,entity); +// if ( nbcells==0 ) return; // e.g. domains have 1 common node + +// map , vector > cellmap; +// map local_offset; +// map distant_offset; + +// //definition of the local offsets for the types present on local +// //and distant domains +// // for a mesh containing 2 triangles and 3 quads +// //local_offset[TRIA3]=0 +// //local_offset[QUAD4]=2 + +// int nb_types_local=(_collection->getMesh())[idomain]-> getNumberOfTypes(entity); +// const MED_EN::medGeometryElement* local_types = (_collection->getMesh())[idomain]->getTypes(entity); +// const int* local_gni = (_collection->getMesh())[idomain]-> getGlobalNumberingIndex(entity); +// for (int i=0; i< nb_types_local; i++) +// { +// local_offset[local_types[i]]=local_gni[i]-1; +// } + +// int nb_types_distant=(_collection->getMesh())[idistant]-> getNumberOfTypes(entity); +// const MED_EN::medGeometryElement* distant_types = (_collection->getMesh())[idistant]->getTypes(entity); +// const int* distant_gni = (_collection->getMesh())[idistant]-> getGlobalNumberingIndex(entity); +// for (int i=0; i< nb_types_distant; i++) +// { +// distant_offset[distant_types[i]]=distant_gni[i]-1; +// } + +// if (nb_types_local==1 && nb_types_distant==1) +// { +// MED_EN::medGeometryElement local_type = (_collection->getMesh())[idomain]->getElementType(entity,1); +// MED_EN::medGeometryElement distant_type = (_collection->getMesh())[idistant]->getElementType(entity,1); +// vector corresp; +// for (int i=0; i key + +// for (int i=0; igetMesh())[idomain]->getElementType(entity,i+1); +// MED_EN::medGeometryElement distant_type = (_collection->getMesh())[idistant]->getElementType(entity,value[icol]); + +// cellmap[make_pair(local_type, distant_type)].push_back(i+1-local_offset[local_type]); +// cellmap[make_pair(local_type, distant_type)].push_back(value[icol]-distant_offset[distant_type]); + +// } +// map , vector >::const_iterator iter; + +// //going through all the (geom,geom) pairs and writing the joints +// for (iter= cellmap.begin(); iter != cellmap.end(); iter++) +// { +// int size= iter->second.size(); +// int *corresp = new int[size]; +// for (int ind=0; ind < size; ind++) +// corresp[ind]=(iter->second)[ind]; +// med_geometrie_element local_geo_elem=(med_geometrie_element)iter->first.first; +// med_geometrie_element distant_geo_elem=(med_geometrie_element)iter->first.second; +// int size_joint=size/2; +// //med_err error = +// MEDjointEcr(fid, mesh_name, joint_name, corresp, size_joint, MED_MAILLE, +// local_geo_elem,MED_MAILLE, distant_geo_elem); +// delete[] corresp; +// } +// // MED v 2.3.1 returns an error code when +// // writing a joint that is already present in the file. +// // Also, it returns an error code if a joint +// // concerns 3D elements. +// // Until these two items are not +// // changed, the following line must be commented out + +// //if (error==-1) throw MEDEXCEPTION("Error filling joint"); + + +// } +// } + +// void MESHCollectionDriver::jointSort(int* elems, int nbelems, bool is_first) +// { +// //filling an ordered structure with the elem ids +// map nodemap; +// if (is_first) +// for (int i=0; i::const_iterator iter=nodemap.begin(); iter!=nodemap.end(); iter++) +// { +// if (is_first) +// { +// *ptr_elems++=iter->first; +// *ptr_elems++=iter->second; +// } +// else +// { +// *ptr_elems++=iter->second; +// *ptr_elems++=iter->first; +// } +// } + + +//} diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.hxx b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.hxx new file mode 100644 index 000000000..e53765d1a --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionDriver.hxx @@ -0,0 +1,86 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MEDPARTITIONER_MESHCOLLECTIONDRIVER_HXX_ +#define MEDPARTITIONER_MESHCOLLECTIONDRIVER_HXX_ + +#include "MEDPARTITIONER.hxx" + +namespace MEDPARTITIONER +{ + class MESHCollection; + class ParaDomainSelector; + + class MEDPARTITIONER_EXPORT MESHCollectionDriver + { + public: + + MESHCollectionDriver(MESHCollection*); + virtual ~MESHCollectionDriver(){} + + virtual int read(char*, ParaDomainSelector* sel=0)=0; + int readSeq(char*,char*); + + virtual void write(char*, ParaDomainSelector* sel=0)=0; +// virtual void readFields (vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) =0; +// virtual void readFields (vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) =0; +// virtual void writeFields(vector *>& filenames, char* fieldname)=0; +// virtual void writeFields(vector *>& filenames, char* fieldname)=0; + +// void readFileStruct(vector & field_names,vector& iternumber,vector & ordernumber,vector & types); + +// int getFieldType(const std::string& fieldname); +// // void exportFamily(vector,MED_EN::medEntityMesh, const string& name); + +// void readLoc2GlobCellConnect(int idomain, const set& loc_domains, ParaDomainSelector* ds, +// vector & loc2glob_corr); + +// int readMeshDimension() const; + + protected: + + void readSubdomain(vector& cellglobal, + vector& faceglobal, + vector& nodeglobal, int idomain); + void writeSubdomain(int idomain,int nbdomain, char*filename, + ParaDomainSelector* domain_selector); + +// void writeElementJoint(medEntityMesh entity , +// int icz, +// int idomain, +// int idistant, +// char* mesh_name, +// char* joint_name, +// med_2_3::med_idt fid ); +// void jointSort(int* elems, int nbelems, bool is_first); + + + + MESHCollection* _collection; + + std::vector _filename; + std::vector _meshname; + + }; + +} + + +#endif /*MESHCOLLECTIONDRIVER_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.H b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.H new file mode 100644 index 000000000..c253b586f --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.H @@ -0,0 +1,66 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MEDSPLITTER_MESHCOLLECTIONMEDASCIIDRIVER_H +#define MEDSPLITTER_MESHCOLLECTIONMEDASCIIDRIVER_H + + +#include +#include +#include +#include + +/*!reads a distributed field + * + * \param fields vector of fields (one field per subdomain) + * \param fieldname name of the field + * \param itnumber number of iteration + * \param ordernumber internal number inside the iteration + * */ +template +void MESHCollectionMedAsciiDriver::_readFields(vector* >& fields,char* fieldname, int itnumber, int ordernumber) +{ + for (int i=0; i<_collection->getMesh().size(); i++) + { + char filename[256]; + strcpy(filename,_filename[i].c_str()); + cout << "maillage : " << filename << " champ : " << fieldname << endl; + // MEDMEM::FIELD* field = new MEDMEM::FIELD(MEDMEM::MED_DRIVER,filename,fieldname,itnumber,ordernumber); + fields.push_back (new MEDMEM::FIELD(MEDMEM::MED_DRIVER,filename,fieldname,itnumber,ordernumber)); + } +} + + +/*!writes a distributed field + * + * \param fields vector of fields (one field per subdomain) + * \param fieldname name of the field + * */ +template +void MESHCollectionMedAsciiDriver::_writeFields(vector* >& fields,char* fieldname) +{ + for (int i=0; i<_collection->getMesh().size(); i++) + { + char filename[256]; + strcpy(filename,_filename[i].c_str()); + int driverid = fields[i]->addDriver(MEDMEM::MED_DRIVER, filename, fieldname); + fields[i]->write(driverid); + } +} + +#endif diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.hxx b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.hxx new file mode 100644 index 000000000..ba45f22c9 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedAsciiDriver.hxx @@ -0,0 +1,76 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MESHCOLLECTIONMEDASCIIDRIVER_HXX_ +#define MESHCOLLECTIONMEDASCIIDRIVER_HXX_ + +#include "MEDPARTITIONER_MESHCollectionDriver.hxx" + +namespace MEDPARTITIONER +{ + class MESHCollection; + + class MESHCollectionMedAsciiDriver:public MESHCollectionDriver + { + public: + + MESHCollectionMedAsciiDriver(MESHCollection*); + virtual ~MESHCollectionMedAsciiDriver(){} + + int read(char*, ParaDomainSelector* sel=0); + + void write(char*, ParaDomainSelector* sel=0); + +// void readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) +// { +// _readFields(filenames,fieldname,itnumber,ordernumber); +// } +// void readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) +// { +// _readFields(filenames,fieldname,itnumber,ordernumber); +// } + +// void writeFields(vector *>& filenames, char* fieldname) +// { +// _writeFields( filenames, fieldname); +// } + +// void writeFields(vector *>& filenames, char* fieldname) +// { +// _writeFields( filenames, fieldname); +// } + + + private : +// template void _readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber); + +// template +// void _writeFields(vector *>& filenames, char* fieldname); + + + + std::string _master_filename; + }; + +} + + +#endif /*MESHCOLLECTIONDRIVER_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.H b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.H new file mode 100644 index 000000000..5d8710b67 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.H @@ -0,0 +1,122 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MEDSPLITTER_MESHCOLLECTIONMEDXMLDRIVER_H +#define MEDSPLITTER_MESHCOLLECTIONMEDXMLDRIVER_H + + +#include +#include +#include +#include + + + +/*!reads a distributed field + * + * \param fields vector of fields (one field per subdomain) + * \param fieldname name of the field + * \param itnumber number of iteration + * \param ordernumber internal number inside the iteration + * */ +template +void MESHCollectionMedXMLDriver::_readFields(vector* >& fields,char* fieldname, int itnumber, int ordernumber) +{ + for (int i=0; i<_collection->getMesh().size(); i++) + { + char filename[256]; + strcpy(filename,_filename[i].c_str()); + cout << "maillage : " << filename << " champ : " << fieldname << endl; + // MEDMEM::FIELD* field = new MEDMEM::FIELD(MEDMEM::MED_DRIVER,filename,fieldname,itnumber,ordernumber); + fields.push_back (new MEDMEM::FIELD(MEDMEM::MED_DRIVER,filename,fieldname,itnumber,ordernumber)); + } +} + + +/*!writes a distributed field + * + * \param fields vector of fields (one field per subdomain) + * \param fieldname name of the field + * */ +template +void MESHCollectionMedXMLDriver::_writeFields(vector* >& fields,char* fieldname) +{ + xmlDocPtr master_doc=xmlParseFile(_master_filename.c_str()); + + if (!master_doc) + throw MEDEXCEPTION("MEDSPLITTER writeFields - Master File does not exist"); + + //number of domains + + xmlXPathContextPtr xpathCtx = xmlXPathNewContext(master_doc); + xmlXPathObjectPtr xpathObj = xmlXPathEvalExpression(BAD_CAST "//mapping/mesh", xpathCtx); + //assuming there is only one mesh in the XML file + xmlNodePtr mesh_node= xpathObj->nodesetval->nodeTab[0]; + + //adds the field to the master file if necessary + bool exist_field =false; + xpathObj = xmlXPathEvalExpression(BAD_CAST "//mapping/mesh/field", xpathCtx); + //assuming there is only one mesh in the XML file + int field_nr = xpathObj->nodesetval->nodeNr; + for (int i=0; inodesetval->nodeTab[i]->properties->children->content, fieldname)==0) + exist_field = true; + } + + xmlNodePtr field_node; + if (!exist_field) + { + field_node = xmlNewChild(mesh_node, 0, BAD_CAST "field",0); + xmlNewProp(field_node,BAD_CAST "name",BAD_CAST fieldname); + } + + + for (int i=0; i<_collection->getMesh().size(); i++) + { + char filename[256]; + strcpy(filename,_filename[i].c_str()); + int driverid = fields[i]->addDriver(MEDMEM::MED_DRIVER, filename, fieldname); + fields[i]->write(driverid); + + //adds the partition to the master file if the field had not been + //added already + if (!exist_field) + { + xmlNodePtr chunk_node= xmlNewChild(field_node,0,BAD_CAST "chunk",0); + char id[8]; + sprintf(id,"%d",i+1); + xmlNewProp(chunk_node,BAD_CAST "subdomain",BAD_CAST id); + //xmlNewProp(chunk_node,BAD_CAST "name", BAD_CAST fieldname); + //xmlNodePtr fieldname_node= + xmlNewChild(chunk_node,0,BAD_CAST "name", BAD_CAST fieldname); + } + } + //Writing file + xmlKeepBlanksDefault(0); + xmlSaveFormatFileEnc(_master_filename.c_str(), master_doc, "UTF-8", 1); + + //Clean up + xmlXPathFreeContext(xpathCtx); + xmlFreeDoc(master_doc); + //xmlCleanupParser(); + +} + +#endif diff --git a/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.hxx b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.hxx new file mode 100644 index 000000000..a36b5623f --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_MESHCollectionMedXMLDriver.hxx @@ -0,0 +1,76 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MESHCOLLECTIONMEDXMLDRIVER_HXX_ +#define MESHCOLLECTIONMEDXMLDRIVER_HXX_ + +#include "MEDPARTITIONER_MESHCollectionDriver.hxx" + +namespace MEDPARTITIONER +{ + class MESHCollection; + + class MESHCollectionMedXMLDriver:public MESHCollectionDriver + { + public: + + MESHCollectionMedXMLDriver(MESHCollection*); + virtual ~MESHCollectionMedXMLDriver(){} + + + int read(char*, ParaDomainSelector* sel=0); + + void write(char*, ParaDomainSelector* sel=0); + +// void readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) +// { +// _readFields(filenames,fieldname,itnumber,ordernumber); +// } +// void readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber) +// { +// _readFields(filenames,fieldname,itnumber,ordernumber); +// } + +// void writeFields(vector *>& filenames, char* fieldname) +// { +// _writeFields( filenames, fieldname); +// } +// void writeFields(vector *>& filenames, char* fieldname) +// { +// _writeFields( filenames, fieldname); +// } + + + + private : + + // template void _readFields(vector *>& filenames, char* fieldname, +// int itnumber, int ordernumber); + +// template +// void _writeFields(vector *>& filenames, char* fieldname); + + std::string _master_filename; + }; + +} + + +#endif /*MESHCOLLECTIONDRIVER_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_METISGraph.cxx b/src/MEDPartitioner/MEDPARTITIONER_METISGraph.cxx new file mode 100644 index 000000000..2c86aacb2 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_METISGraph.cxx @@ -0,0 +1,120 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +#ifdef ENABLE_PARMETIS +#include +#endif +extern "C" { +#include +} +#include "MEDPARTITIONER_METISGraph.hxx" +#include "MEDPARTITIONER_ParaDomainSelector.hxx" +#include "MEDMEM_Exception.hxx" +using namespace MEDPARTITIONER; + +METISGraph::METISGraph():Graph() +{ +} + +METISGraph::METISGraph(MEDPARTITIONER::MEDSKYLINEARRAY* graph, int* edgeweight) + :Graph(graph,edgeweight) +{ +} + +METISGraph::~METISGraph() +{ +} + +void METISGraph::partGraph(int ndomain, + const std::string& options_string, + ParaDomainSelector* parallelizer) +{ + // number of graph vertices + int n = m_graph->getNumberOf(); + + //graph + int * xadj=const_cast(m_graph->getIndex()); + int * adjncy = const_cast(m_graph->getValue()); + //constraints + int * vwgt=m_cellweight; + int * adjwgt=m_edgeweight; + int wgtflag=(m_edgeweight!=0)?1:0+(m_cellweight!=0)?2:0; + + //base 0 or 1 + int base=0; + + //ndomain + int nparts = ndomain; + + //options + int options[5]={0,0,0,0,0}; + + // output parameters + int edgecut; + int* partition = new int[n]; + + if (nparts >1) + { + if ( parallelizer ) + { +#ifdef ENABLE_PARMETIS + // distribution of vertices of the graph among the processors + int * vtxdist = parallelizer ? parallelizer->getNbVertOfProcs() : 0; + MPI_Comm comm = MPI_COMM_WORLD; + + ParMETIS_PartKway( vtxdist, xadj, adjncy, vwgt, adjwgt, &wgtflag, + &base, &nparts, options, &edgecut, partition, &comm ); +#else + throw MEDMEM::MEDEXCEPTION("ParMETIS is not available. Check your products, please."); +#endif + } + else + { + + if (options_string != "k") + METIS_PartGraphRecursive(&n, xadj, adjncy, vwgt, adjwgt, &wgtflag, + &base, &nparts, options, &edgecut, partition); + else + METIS_PartGraphKway(&n, xadj, adjncy, vwgt, adjwgt, &wgtflag, + &base, &nparts, options, &edgecut, partition); + } + } + else + { + for (int i=0; i index(n+1); + std::vector value(n); + index[0]=0; + for (int i=0; i +namespace MEDPARTITIONER { + class MEDSKYLINEARRAY; + class MEDPARTITIONER_EXPORT METISGraph:public Graph + { + public: + METISGraph(); + METISGraph(MEDPARTITIONER::MEDSKYLINEARRAY*, int* edgeweight=0); + virtual ~METISGraph(); + void partGraph(int ndomain, const std::string& options_string="", ParaDomainSelector* sel=0); + //private: + // const MEDMEM::MEDSKYLINEARRAY* m_graph; + }; +} +#endif /*METISGRAPH_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.cxx b/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.cxx new file mode 100644 index 000000000..1c0a660c5 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.cxx @@ -0,0 +1,667 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// File : MEDSPLITTER_ParaDomainSelector.cxx +// Created : Wed Jun 24 12:39:59 2009 +// Author : Edward AGAPOV (eap) + +#include "MEDCouplingUMesh.hxx" +#include "MEDMEM_Exception.hxx" +#include "MEDPARTITIONER_ParaDomainSelector.hxx" + +#include "MEDPARTITIONER_UserGraph.hxx" +#include "MEDPARTITIONER_JointExchangeData.hxx" + +//#include +#include + +#include +#include + +#ifdef HAVE_MPI2 +#include +#endif + +#ifndef WIN32 +#include +#endif + +using namespace MEDPARTITIONER; +using namespace MED_EN; +using namespace std; + +//================================================================================ +/*! + * \brief Constructor. Find out my rank and world size + */ +//================================================================================ + +ParaDomainSelector::ParaDomainSelector(bool mesure_memory) + :_rank(0),_world_size(1), _nb_result_domains(-1), _mesure_memory(mesure_memory), + _init_time(0.0), _init_memory(0), _max_memory(0) +{ +#ifdef HAVE_MPI2 + MPI_Comm_size(MPI_COMM_WORLD,&_world_size) ; + MPI_Comm_rank(MPI_COMM_WORLD,&_rank) ; + _init_time = MPI_Wtime(); +#endif + evaluateMemory(); +} + +//================================================================================ +/*! + * \brief Destructor + */ +//================================================================================ + +ParaDomainSelector::~ParaDomainSelector() +{ +} + +//================================================================================ +/*! + * \brief Return true if is running on different hosts + */ +//================================================================================ + +bool ParaDomainSelector::isOnDifferentHosts() const +{ + evaluateMemory(); + if ( _world_size < 2 ) + return false; + +#ifdef HAVE_MPI2 + char name_here[ MPI_MAX_PROCESSOR_NAME+1 ], name_there[ MPI_MAX_PROCESSOR_NAME+1 ]; + int size; + MPI_Get_processor_name( name_here, &size); + + int next_proc = (rank() + 1) % nbProcs(); + int prev_proc = (rank() - 1 + nbProcs()) % nbProcs(); + int tag = 1111111; + + MPI_Status status; + MPI_Sendrecv((void*)&name_here[0], MPI_MAX_PROCESSOR_NAME, MPI_CHAR, next_proc, tag, + (void*)&name_there[0], MPI_MAX_PROCESSOR_NAME, MPI_CHAR, prev_proc, tag, + MPI_COMM_WORLD, &status); + return string(name_here) != string(name_there); +#endif +} + +//================================================================================ +/*! + * \brief Return true if the domain with domainIndex is to be loaded on this proc + * \param domainIndex - index of mesh domain + * \retval bool - to load or not + */ +//================================================================================ + +bool ParaDomainSelector::isMyDomain(int domainIndex) const +{ + evaluateMemory(); + return (_rank == getProccessorID( domainIndex )); +} + +//================================================================================ +/*! + * \brief Return processor id where the domain with domainIndex resides + * \param domainIndex - index of mesh domain + * \retval int - processor id + */ +//================================================================================ + +int ParaDomainSelector::getProccessorID(int domainIndex) const +{ + evaluateMemory(); + return ( domainIndex % _world_size ); +} + +//================================================================================ +/*! + * \brief Gather info on nb of entities on each processor and return total nb. + * + * Is called + * 1) for MED_CELL to know global id shift for domains at graph construction; + * 2) for sub-entity to know total nb of sub-entities before creating those of joints + */ +//================================================================================ + +int ParaDomainSelector::gatherNbOf( + //MED_EN::medEntityMesh entity, + const vector& domain_meshes) +{ + evaluateMemory(); + + // get nb of elems of each domain mesh + int nb_domains = domain_meshes.size(); + vector nb_elems( nb_domains, 0 ); + for ( int i = 0; i < nb_domains; ++i ) + if ( domain_meshes[i] ) + nb_elems[i] = domain_meshes[i]->getNumberOfCells(); + + // receive nb of elems from other procs + vector all_nb_elems( nb_domains ); +#ifdef HAVE_MPI2 + MPI_Allreduce((void*)&nb_elems[0], (void*)&all_nb_elems[0], nb_domains, + MPI_INT, MPI_SUM, MPI_COMM_WORLD); +#endif + int total_nb = std::accumulate( all_nb_elems.begin(), all_nb_elems.end(), 0 ); + + vector& elem_shift_by_domain = _cell_shift_by_domain; + + // fill elem_shift_by_domain + + vector< int > ordered_nbs, domain_order( nb_domains ); + ordered_nbs.push_back(0); + for ( int iproc = 0; iproc < nbProcs(); ++iproc ) + for ( int idomain = 0; idomain < nb_domains; ++idomain ) + if ( getProccessorID( idomain ) == iproc ) + { + domain_order[ idomain ] = ordered_nbs.size() - 1; + ordered_nbs.push_back( ordered_nbs.back() + all_nb_elems[idomain] ); + } + elem_shift_by_domain.resize( nb_domains+1, 0 ); + for ( int idomain = 0; idomain < nb_domains; ++idomain ) + elem_shift_by_domain[ idomain ] = ordered_nbs[ domain_order[ idomain ]]; + + elem_shift_by_domain.back() = ordered_nbs.back(); // to know total nb of elements + + // if ( entity == MED_CELL ) + { + // fill _nb_vert_of_procs + _nb_vert_of_procs.resize( _world_size+1, 0 ); + for ( int i = 0; i < nb_domains; ++i ) + { + int rank = getProccessorID( i ); + _nb_vert_of_procs[ rank+1 ] += all_nb_elems[ i ]; + } + _nb_vert_of_procs[0] = 1; // base = 1 + for ( int i = 1; i < _nb_vert_of_procs.size(); ++i ) + _nb_vert_of_procs[ i ] += _nb_vert_of_procs[ i-1 ]; // to CSR format + } +// else +// { +// // to compute global ids of faces in joints +// //_total_nb_faces = total_nb; +// } + +// if ( !_rank) { +// MEDMEM::STRING out("_nb_vert_of_procs :"); +// for ( int i = 0; i < _nb_vert_of_procs.size(); ++i ) +// out << _nb_vert_of_procs[i] << " "; +// std::cout << out << std::endl; +// } +// if ( !_rank) { +// MEDMEM::STRING out("elem_shift_by_domain :"); +// for ( int i = 0; i < elem_shift_by_domain.size(); ++i ) +// out << elem_shift_by_domain[i] << " "; +// std::cout << out << std::endl; +// } + + evaluateMemory(); + + return total_nb; +} + +//================================================================================ +/*! + * \brief Return distribution of the graph vertices among the processors + * \retval int* - array conatining nb of vertices on all processors + * + * gatherNbOf( MED_CELL ) must be called before. + * The result array is to be used as the first arg of ParMETIS_V3_PartKway() and + * is freed by ParaDomainSelector. + */ +//================================================================================ + +#define gatherNbOf_NOT_CALLED(meth) throw MED_EXCEPTION \ +("ParaDomainSelector::" #meth "(): gatherNbOf( MED_CELL ) must be called before") + +int* ParaDomainSelector::getNbVertOfProcs() const +{ + evaluateMemory(); + if ( _nb_vert_of_procs.empty() ) + gatherNbOf_NOT_CALLED(getNbVertOfProcs); + + return (int*) & _nb_vert_of_procs[0]; +} +//================================================================================ +/*! + * \brief Return nb of cells in domains with lower index. + * + * gatherNbOf( MED_CELL ) must be called before. + * Result added to local id on given domain gives id in the whole distributed mesh + */ +//================================================================================ + +int ParaDomainSelector::getDomainShift(int domainIndex) const +{ + evaluateMemory(); + if ( _cell_shift_by_domain.empty() ) + gatherNbOf_NOT_CALLED(getDomainShift); + + return _cell_shift_by_domain[ domainIndex ]; +} + +//================================================================================ +/*! + * \brief Return nb of cells on processors with lower rank. + * + * gatherNbOf( MED_CELL ) must be called before. + * Result added to global id on this processor gives id in the whole distributed mesh + */ +//================================================================================ + +int ParaDomainSelector::getProcShift() const +{ + evaluateMemory(); + if ( _nb_vert_of_procs.empty() ) + gatherNbOf_NOT_CALLED(getProcShift); + + return _nb_vert_of_procs[_rank]-1; +} + +//================================================================================ +/*! + * \brief Gather graphs from all processors into one + */ +//================================================================================ + +auto_ptr ParaDomainSelector::gatherGraph(const Graph* graph) const +{ + Graph* glob_graph = 0; + + evaluateMemory(); +#ifdef HAVE_MPI2 + + // --------------- + // Gather indices + // --------------- + + vector index_size_of_proc( nbProcs() ); // index sizes - 1 + for ( int i = 1; i < _nb_vert_of_procs.size(); ++i ) + index_size_of_proc[i-1] = _nb_vert_of_procs[ i ] - _nb_vert_of_procs[ i-1 ]; + + int index_size = 1 + _cell_shift_by_domain.back(); + int* graph_index = new int[ index_size ]; + const int* index = graph->getGraph()->getIndex(); + int* proc_index_displacement = (int*) & _nb_vert_of_procs[0]; + + MPI_Allgatherv((void*) (index+1), // send local index except first 0 (or 1) + index_size_of_proc[_rank], // index size on this proc + MPI_INT, + (void*) graph_index, // receive indices + & index_size_of_proc[0], // index size on each proc + proc_index_displacement, // displacement of each proc index + MPI_INT, + MPI_COMM_WORLD); + graph_index[0] = index[0]; // it is not overwritten thanks to proc_index_displacement[0]==1 + + // get sizes of graph values on each proc by the got indices of graphs + vector< int > value_size_of_proc( nbProcs() ), proc_value_displacement(1,0); + for ( int i = 0; i < nbProcs(); ++i ) + { + if ( index_size_of_proc[i] > 0 ) + value_size_of_proc[i] = graph_index[ proc_index_displacement[ i+1 ]-1 ] - graph_index[0]; + else + value_size_of_proc[i] = 0; + proc_value_displacement.push_back( proc_value_displacement.back() + value_size_of_proc[i] ); + } + + // update graph_index + for ( int i = 1; i < nbProcs(); ++i ) + { + int shift = graph_index[ proc_index_displacement[i]-1 ]-graph_index[0]; + for ( int j = proc_index_displacement[i]; j < proc_index_displacement[i+1]; ++j ) + graph_index[ j ] += shift; + } + + // -------------- + // Gather values + // -------------- + + int value_size = graph_index[ index_size-1 ] - graph_index[ 0 ]; + int* graph_value = new int[ value_size ]; + const int* value = graph->getGraph()->getValue(); + + MPI_Allgatherv((void*) value, // send local value + value_size_of_proc[_rank], // value size on this proc + MPI_INT, + (void*) graph_value, // receive values + & value_size_of_proc[0], // value size on each proc + & proc_value_displacement[0], // displacement of each proc value + MPI_INT, + MPI_COMM_WORLD); + + // ----------------- + // Gather partition + // ----------------- + + int * partition = new int[ _cell_shift_by_domain.back() ]; + const int* part = graph->getPart(); + + MPI_Allgatherv((void*) part, // send local partition + index_size_of_proc[_rank], // index size on this proc + MPI_INT, + (void*)(partition-1), // -1 compensates proc_index_displacement[0]==1 + & index_size_of_proc[0], // index size on each proc + proc_index_displacement, // displacement of each proc index + MPI_INT, + MPI_COMM_WORLD); + + // ----------- + // Make graph + // ----------- + + MEDMEM::MEDSKYLINEARRAY* array = + new MEDMEM::MEDSKYLINEARRAY( index_size-1, value_size, graph_index, graph_value, true ); + + glob_graph = new UserGraph( array, partition, index_size-1 ); + + evaluateMemory(); + + delete [] partition; + +#endif // HAVE_MPI2 + + return auto_ptr( glob_graph ); +} + +//================================================================================ +/*! + * \brief Sets global numbering for the entity. + * + * This method must be once called for MED_CELL before exchangeJoint() calls + */ +//================================================================================ + +// void ParaDomainSelector::gatherEntityTypesInfo(vector& domain_meshes, +// MED_EN::medEntityMesh entity) +// { +// const list & all_types = meshEntities[ entity ]; + +// evaluateMemory(); + +// // Put nb of elements of the entity of all domains in one vector +// // and by the way find out mesh dimension + +// vector nb_of_type( domain_meshes.size() * all_types.size(), 0 ); +// int mesh_dim = -1, space_dim = -1; + +// for ( int idomain = 0; idomain < domain_meshes.size(); ++idomain ) +// { +// if ( !isMyDomain(idomain)) continue; + +// int* domain_nbs = & nb_of_type[ idomain * all_types.size()]; + +// list::const_iterator type = all_types.begin(); +// for ( int t = 0; type != all_types.end(); ++t,++type ) +// domain_nbs[t] = domain_meshes[idomain]->getNumberOfCells(); + +// int i_mesh_dim = domain_meshes[idomain]->getMeshDimension(); +// int i_space_dim = domain_meshes[idomain]->getSpaceDimension(); +// if ( mesh_dim < i_mesh_dim && i_mesh_dim <= 3 ) +// mesh_dim = i_mesh_dim; +// if ( space_dim < i_space_dim && i_space_dim <= 3 ) +// space_dim = i_space_dim; +// } + +// // Receive nbs from other procs + +// vector< int > nb_recv( nb_of_type.size() ); +// #ifdef HAVE_MPI2 +// MPI_Allreduce((void*)&nb_of_type[0], (void*)&nb_recv[0], nb_of_type.size(), +// MPI_INT, MPI_SUM, MPI_COMM_WORLD); +// #endif + +// // Set info to meshes of distant domains + +// for ( int idomain = 0; idomain < domain_meshes.size(); ++idomain ) +// { +// if ( isMyDomain(idomain)) continue; + +// MEDMEM::MESHING* meshing = (MEDMEM::MESHING*) domain_meshes[idomain]; +// if ( meshing->getMeshDimension() < mesh_dim ) +// { +// meshing->setMeshDimension( mesh_dim ); +// meshing->setSpaceDimension( space_dim ); +// } + +// vector< medGeometryElement > types; +// vector< int > nb_elems; + +// int* domain_nbs = & nb_recv[ idomain * all_types.size()]; + +// list::const_iterator type = all_types.begin(); +// for ( int t = 0; type != all_types.end(); ++t,++type ) +// { +// if ( domain_nbs[t]==0 ) continue; +// types.push_back( *type ); +// nb_elems.push_back( domain_nbs[t] ); +// } +// meshing->setNumberOfTypes( types.size(), entity ); +// if ( !types.empty() ) +// { +// meshing->setTypes ( & types[0], entity ); +// meshing->setNumberOfElements( & nb_elems[0], entity ); +// } +// } +// evaluateMemory(); +// } + +//================================================================================ +/*! + * \brief Set nb of cell/cell pairs in a joint between domains + */ +//================================================================================ + +void ParaDomainSelector::setNbCellPairs( int nb_cell_pairs, int dist_domain, int loc_domain ) +{ + // This method is needed for further computing global numbers of faces in joint. + // Store if both domains are on this proc else on one of procs only + if ( isMyDomain( dist_domain ) || dist_domain < loc_domain ) + { + if ( _nb_cell_pairs_by_joint.empty() ) + _nb_cell_pairs_by_joint.resize( _nb_result_domains*(_nb_result_domains+1), 0); + + int joint_id = jointId( loc_domain, dist_domain ); + _nb_cell_pairs_by_joint[ joint_id ] = nb_cell_pairs; + } + evaluateMemory(); +} + +//================================================================================ +/*! + * \brief Return nb of cell/cell pairs in a joint between domains on different procs + */ +//================================================================================ + +int ParaDomainSelector::getNbCellPairs( int dist_domain, int loc_domain ) const +{ + evaluateMemory(); + + int joint_id = jointId( loc_domain, dist_domain ); + return _nb_cell_pairs_by_joint[ joint_id ]; +} + +//================================================================================ +/*! + * \brief Gather size of each joint + */ +//================================================================================ + +void ParaDomainSelector::gatherNbCellPairs() +{ + const char* LOC = "MEDPARTITIONER::ParaDomainSelector::gatherNbCellPairs(): "; + if ( _nb_cell_pairs_by_joint.empty() ) + _nb_cell_pairs_by_joint.resize( _nb_result_domains*(_nb_result_domains+1), 0); + evaluateMemory(); + + vector< int > send_buf = _nb_cell_pairs_by_joint; +#ifdef HAVE_MPI2 + MPI_Allreduce((void*)&send_buf[0], + (void*)&_nb_cell_pairs_by_joint[0], + _nb_cell_pairs_by_joint.size(), + MPI_INT, MPI_SUM, MPI_COMM_WORLD); +#endif + // check that the set nbs of cell pairs are correct, + // namely that each joint is treated on one proc only + for ( int j = 0; j < _nb_cell_pairs_by_joint.size(); ++j ) + if ( _nb_cell_pairs_by_joint[j] != send_buf[j] && send_buf[j]>0 ) + throw MED_EXCEPTION("invalid nb of cell pairs"); +} + +//================================================================================ +/*! + * \brief Send-receive joint data + */ +//================================================================================ + +// void ParaDomainSelector::exchangeJoint( JointExchangeData* joint ) const +// { +// vector send_data, recv_data( joint->serialize( send_data )); + +// int dest = getProccessorID( joint->distantDomain() ); +// int tag = 1001 + jointId( joint->localDomain(), joint->distantDomain() ); + +// #ifdef HAVE_MPI2 +// MPI_Status status; +// MPI_Sendrecv((void*)&send_data[0], send_data.size(), MPI_INT, dest, tag, +// (void*)&recv_data[0], recv_data.size(), MPI_INT, dest, tag, +// MPI_COMM_WORLD, &status); +// #endif + +// joint->deserialize( recv_data ); +// } + +//================================================================================ +/*! + * \brief Return the first global id of sub-entity for the joint + */ +//================================================================================ + +int ParaDomainSelector::getFisrtGlobalIdOfSubentity( int loc_domain, int dist_domain ) const +{ + // total_nb_faces includes faces existing before creation of joint faces + // (got in gatherNbOf( MED_FACE )). + evaluateMemory(); + + int total_nb_faces = _face_shift_by_domain.empty() ? 0 : _face_shift_by_domain.back(); + int id = total_nb_faces + 1; + + if ( _nb_cell_pairs_by_joint.empty() ) + throw MED_EXCEPTION("MEDPARTITIONER::ParaDomainSelector::getFisrtGlobalIdOfSubentity(), " + "gatherNbCellPairs() must be called before"); + int joint_id = jointId( loc_domain, dist_domain ); + for ( int j = 0; j < joint_id; ++j ) + id += _nb_cell_pairs_by_joint[ j ]; + + return id; +} + +//================================================================================ +/*! + * \brief Send-receive local ids of joint faces + */ +//================================================================================ + +int* ParaDomainSelector::exchangeSubentityIds( int loc_domain, int dist_domain, + const vector& loc_ids_here ) const +{ + int* loc_ids_dist = new int[ loc_ids_here.size()]; + int dest = getProccessorID( dist_domain ); + int tag = 2002 + jointId( loc_domain, dist_domain ); +#ifdef HAVE_MPI2 + MPI_Status status; + MPI_Sendrecv((void*)&loc_ids_here[0], loc_ids_here.size(), MPI_INT, dest, tag, + (void*) loc_ids_dist, loc_ids_here.size(), MPI_INT, dest, tag, + MPI_COMM_WORLD, &status); +#endif + evaluateMemory(); + + return loc_ids_dist; +} + +//================================================================================ +/*! + * \brief Return identifier for a joint + */ +//================================================================================ + +int ParaDomainSelector::jointId( int local_domain, int distant_domain ) const +{ + evaluateMemory(); + if (_nb_result_domains < 0) + throw MED_EXCEPTION("ParaDomainSelector::jointId(): setNbDomains() must be called before()"); + + if ( local_domain < distant_domain ) + swap( local_domain, distant_domain ); + return local_domain * _nb_result_domains + distant_domain; +} + +//================================================================================ +/*! + * \brief Return domain order so that first go domains on proc 0 and so n + */ +//================================================================================ + +// int ParaDomainSelector::getDomianOrder(int idomain, int nb_domains) const +// { +// return nb_domains / nbProcs() * getProccessorID( idomain ) + idomain / nbProcs(); +// } + +//================================================================================ +/*! + * \brief Return time passed from construction in seconds + */ +//================================================================================ + +double ParaDomainSelector::getPassedTime() const +{ +#ifdef HAVE_MPI2 + return MPI_Wtime() - _init_time; +#else + return 0.0; +#endif +} + +//================================================================================ +/*! + * \brief Evaluate current memory usage and return the maximal one in KB + */ +//================================================================================ + +int ParaDomainSelector::evaluateMemory() const +{ + if ( _mesure_memory ) + { + int used_memory = 0; +#ifndef WIN32 + struct sysinfo si; + int err = sysinfo( &si ); + if ( !err ) + used_memory = + (( si.totalram - si.freeram + si.totalswap - si.freeswap ) * si.mem_unit ) / 1024; +#endif + if ( used_memory > _max_memory ) + ((ParaDomainSelector*) this)->_max_memory = used_memory; + + if ( !_init_memory ) + ((ParaDomainSelector*) this)->_init_memory = used_memory; + } + return _max_memory - _init_memory; +} diff --git a/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.hxx b/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.hxx new file mode 100644 index 000000000..ee502b4b2 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_ParaDomainSelector.hxx @@ -0,0 +1,151 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// File : MEDPARTITIONER_ParaDomainSelector.hxx +// Created : Wed Jun 24 12:18:58 2009 +// Author : Edward AGAPOV (eap) + + +#ifndef __MEDPARTITIONER_ParaDomainSelector_HXX__ +#define __MEDPARTITIONER_ParaDomainSelector_HXX__ + +#include "MEDPARTITIONER.hxx" + +#include + +#include +#include + +namespace ParaMEDMEM +{ + class MEDCouplingUMesh; +} + + +namespace MEDPARTITIONER +{ + class Graph; + class JointExchangeData; + +/*! + * \brief Communication helper in parallel mode + */ +class MEDPARTITIONER_EXPORT ParaDomainSelector +{ +public: + + ParaDomainSelector(bool mesure_memory=false); + ~ParaDomainSelector(); + + //!< return processor rank + int rank() const { return _rank; } + + //!< return number of processors + int nbProcs() const { return _world_size; } + + // Return true if is running on different hosts + bool isOnDifferentHosts() const; + + // Return true if the domain with domainIndex is to be loaded on this proc + bool isMyDomain(int domainIndex) const; + + // Return processor id where the domain with domainIndex resides + int getProccessorID(int domainIndex) const; + + + //!< Set nb of required domains. (Used to sort joints via jointId()) + void setNbDomains(int nb) { _nb_result_domains = nb; } + + // Return identifier for a joint + int jointId( int local_domain, int distant_domain ) const; + + // Return domain order + //int getDomianOrder(int domainIndex, int nb_domains) const; + + + // Collect nb of entities on procs and return total nb + int gatherNbOf( + //MED_EN::medEntityMesh entity, +const std::vector& domain_meshes); + + // Return distribution of the graph vertices among the processors + int* getNbVertOfProcs() const; + + // Return nb of cells on processors with lower rank + int getProcShift() const; + + // Return nb of cells in domains with lower index + int getDomainShift(int domainIndex) const; + +// // Return nb of sub-entities in domains with lower index +// int getDomainSubentityShift(int domainIndex) const; + + // Gather graphs from all processors into one + std::auto_ptr gatherGraph(const Graph* graph) const; + + // Set types and number of elements of the entity to all meshes +// void gatherEntityTypesInfo(std::vector& domain_meshes, +// MED_EN::medEntityMesh entity); + + // Set nb of cell/cell pairs in a joint between domains + void setNbCellPairs( int nb_cell_pairs, int dist_domain, int loc_domain ); + + // Gather size of each proc/proc joint + void gatherNbCellPairs(); + + // Return nb of cell/cell pairs in a joint between domains on different procs + int getNbCellPairs( int dist_domain, int loc_domain ) const; + + // Send-receive joint data + // void exchangeJoint( JointExchangeData* joint ) const; + + // Return the first global id of sub-entity for the joint + int getFisrtGlobalIdOfSubentity( int loc_domain, int dist_domain ) const; + + // Send-receive local ids of joint faces + int* exchangeSubentityIds( int loc_domain, int dist_domain, + const std::vector& loc_ids_here ) const; + + // Return time passed from construction in seconds + double getPassedTime() const; + + // Evaluate current memory usage and return the maximal one in KB + int evaluateMemory() const; + +private: + + int _rank, _world_size; // my rank and nb of processors + + int _nb_result_domains; // required nb of domains + + //int _total_nb_faces; // nb of faces in the whole mesh without proc/proc joint faces + + std::vector< int > _nb_cell_pairs_by_joint; + + std::vector< int > _nb_vert_of_procs; // graph vertices + std::vector< int > _cell_shift_by_domain; + std::vector< int > _face_shift_by_domain; + + double _init_time; + bool _mesure_memory; + int _init_memory, _max_memory; +}; + +} + +#endif diff --git a/src/MEDPartitioner/MEDPARTITIONER_ParallelTopology.cxx b/src/MEDPartitioner/MEDPARTITIONER_ParallelTopology.cxx new file mode 100644 index 000000000..1735e3a4a --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_ParallelTopology.cxx @@ -0,0 +1,1250 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#include +#include +#include +#ifndef WNT +# include +#else +# include +#endif + +#include "MEDPARTITIONER_MESHCollection.hxx" +#include "MEDPARTITIONER_Topology.hxx" +#include "MEDPARTITIONER_Graph.hxx" +#include "MEDPARTITIONER_ParallelTopology.hxx" +#include "MEDPARTITIONER_ConnectZone.hxx" + +#include "MEDCouplingUMesh.hxx" +#include "MEDMEM_Exception.hxx" +#include "MEDMEM_Utilities.hxx" + +#ifndef WNT +using namespace __gnu_cxx; +#else +using namespace std; +#endif + +using namespace MEDPARTITIONER; + +//empty constructor +ParallelTopology::ParallelTopology():m_nb_domain(0),m_mesh_dimension(0) +{} + +//!constructing topology according to mesh collection +ParallelTopology::ParallelTopology(const vector& meshes, + const vector& cz, + vector& cellglobal, + vector& nodeglobal, + vector& faceglobal):m_nb_domain(meshes.size())/*,m_mesh_dimension(meshes[0]->getMeshDimension())*/ +{ + + int index_global=0; + int index_node_global=0; + int index_face_global=0; + + m_nb_cells.resize(m_nb_domain); + m_nb_nodes.resize(m_nb_domain); + // m_nb_faces.resize(m_nb_domain); + + m_loc_to_glob.resize(m_nb_domain); + m_node_loc_to_glob.resize(m_nb_domain); + // m_face_loc_to_glob.resize(m_nb_domain); + + //MED_EN::medEntityMesh constituent_entity; + + bool parallel_mode = false; + for (int idomain=0; !parallel_mode && idomaingetMeshDimension(); + //constituent_entity = (m_mesh_dimension == 3 ? MED_EN::MED_FACE : MED_EN::MED_EDGE ); + + //creating cell maps + m_nb_cells[idomain]=meshes[idomain]->getNumberOfCells(); + // cout << "Nb cells (domain "< ("< ("<getNumberOfNodes()); + for (int i=0; igetNumberOfNodes(); i++) + { + m_node_glob_to_loc.insert(make_pair(i,make_pair(0,i))); + //m_node_loc_to_glob.insert(make_pair(make_pair(0,i+1), i+1)); + m_node_loc_to_glob[0][i]=i; + } + m_nb_total_nodes=meshes[idomain]->getNumberOfNodes(); + m_nb_nodes[0]=m_nb_total_nodes; + + // meshes[idomain]->getConnectivity(MED_EN::MED_FULL_INTERLACE, MED_EN::MED_DESCENDING, MED_EN::MED_CELL, MED_EN::MED_ALL_ELEMENTS); + // int nbfaces=meshes[idomain]->getNumberOfElementsWithPoly(constituent_entity,MED_EN::MED_ALL_ELEMENTS); +// m_face_loc_to_glob[idomain].resize(nbfaces); +// for (int i=0; igetNumberOfNodes(); + hash_map > local2distant; + m_node_loc_to_glob[idomain].resize(m_nb_nodes[idomain]); + for (int icz=0; iczgetLocalDomainNumber() == idomain && + cz[icz]->getLocalDomainNumber()>cz[icz]->getDistantDomainNumber()) + { + int nb_node= cz[icz]->getNodeNumber(); + const int* node_corresp=cz[icz]->getNodeCorrespValue(); + int distant_ip = cz[icz]->getDistantDomainNumber(); + for (int i=0; i< nb_node; i++) + { + int local= node_corresp[i*2]; + int distant = node_corresp[i*2+1]; + local2distant.insert(make_pair(local, make_pair(distant_ip,distant))); + } + } + } + // setting mappings for all nodes + if (nodeglobal[idomain]==0) + { + for (int inode=0; inodesecond).first; + int distant = (local2distant.find(inode)->second).second; + //int global_number=m_loc_to_glob[make_pair(ip,distant)]; + int global_number=m_loc_to_glob[ip][distant]; + m_node_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,inode))); + //m_node_loc_to_glob[make_pair(idomain,inode+1)]=global_number; + m_node_loc_to_glob[idomain][inode]=global_number; + } + } + } + //using former node numbering + else + {// cout << "("<"<getConnectivity(MED_EN::MED_FULL_INTERLACE, MED_EN::MED_DESCENDING, MED_EN::MED_CELL, MED_EN::MED_ALL_ELEMENTS); +// m_nb_faces[idomain]=meshes[idomain]->getNumberOfElementsWithPoly(constituent_entity,MED_EN::MED_ALL_ELEMENTS); +// MESSAGE_MED ("Nb faces domain " << idomain<getLocalDomainNumber() == idomain && +// cz[icz]->getLocalDomainNumber()>cz[icz]->getDistantDomainNumber()) +// { +// int nb_face= cz[icz]->getFaceNumber(); +// const int* face_corresp=cz[icz]->getFaceCorrespValue(); +// int distant_ip = cz[icz]->getDistantDomainNumber(); +// for (int i=0; i< nb_face; i++) +// { +// int local= face_corresp[i*2]; +// int distant = face_corresp[i*2+1]; +// local2distant.insert(make_pair(local, make_pair(distant_ip,distant))); +// } +// } +// } +// // setting mappings for all faces +// if (faceglobal[idomain]==0) +// { +// for (int iface=0; ifacesecond).first; +// int distant = (local2distant.find(iface+1)->second).second; +// //int global_number=m_loc_to_glob[make_pair(ip,distant)]; +// int global_number=m_loc_to_glob[ip][distant-1]; +// m_face_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,iface+1))); +// //m_face_loc_to_glob[make_pair(idomain,iface+1)]=global_number; +// m_face_loc_to_glob[idomain][iface]=global_number; +// } +// } +// } +// //using former face numbering +// else +// { +// for (int iface=0; ifacenbVertices()), + m_mesh_dimension(mesh_dimension) +{ + m_nb_cells.resize(m_nb_domain); + m_nb_nodes.resize(m_nb_domain); + m_nb_faces.resize(m_nb_domain); + + m_loc_to_glob.resize(m_nb_domain); + m_node_loc_to_glob.resize(m_nb_domain); + m_face_loc_to_glob.resize(m_nb_domain); + + // used in parallel mode only + m_cell_loc_to_glob_fuse.resize(m_nb_domain); + m_face_loc_to_glob_fuse.resize(m_nb_domain); + + for (int i=0; i getPart(); + m_nb_total_cells= graph->nbVertices(); + + for (int icell=0; icell local_node = m_node_glob_to_loc.find(node_list[i])->second; + ip[i]=local_node.first; + local[i]=local_node.second; + } +} + +/*!Converts a list of global node numbers on domain ip + * to a distributed array with local cell numbers. + * + * If a node in the list is represented on several domains, + * only the value with domain ip is returned + * + * */ +void ParallelTopology::convertGlobalNodeList(const int* node_list, int nbnode, int* local, int ip) +{ + if (m_node_glob_to_loc.empty()) + throw MEDMEM::MEDEXCEPTION("convertGlobalNodeList - Node mapping has not yet been built"); + + for (int i=0; i< nbnode; i++) + { + typedef hash_multimap >::iterator mmiter; + pair range=m_node_glob_to_loc.equal_range(node_list[i]); + for (mmiter it=range.first; it !=range.second; it++) + { + int ipfound=(it->second).first; + if (ipfound==ip) + local[i]=(it->second).second; + } + } +} + +/*!Converts a list of global node numbers + * to a distributed array with local cell numbers. + * + * If a node in the list is represented on several domains, + * all the values are put in the array + * */ +void ParallelTopology::convertGlobalNodeListWithTwins(const int* node_list, int nbnode, int*& local, int*& ip,int*& full_array, int& size) +{ + if (m_node_glob_to_loc.empty()) + throw MEDMEM::MEDEXCEPTION("convertGlobalNodeList - Node mapping has not yet been built"); + + size=0; + for (int i=0; i< nbnode; i++) + { + int count= m_node_glob_to_loc.count(node_list[i]); + // if (count > 1) + // cout << "noeud " << node_list[i]<< " doublon d'ordre " << count< >::iterator mmiter; + pair range=m_node_glob_to_loc.equal_range(node_list[i]); + for (mmiter it=range.first; it !=range.second; it++) + { + ip[index]=(it->second).first; + local[index]=(it->second).second; + full_array [index]=node_list[i]; + index++; + } + + } +} + +/*!Converts a list of global face numbers + * to a distributed array with local face numbers. + * + * If a face in the list is represented on several domains, + * all the values are put in the array + * */ +void ParallelTopology::convertGlobalFaceListWithTwins(const int* face_list, int nbface, int*& local, int*& ip, int*& full_array,int& size) +{ + size=0; + for (int i=0; i< nbface; i++) + { + //int count = m_face_glob_to_loc.count(face_list[i]); + //if (count >1) MESSAGE_MED("face en doublon "< >::iterator mmiter; + pair range=m_face_glob_to_loc.equal_range(face_list[i]); + for (mmiter it=range.first; it !=range.second; it++) + { + ip[index]=(it->second).first; + local[index]=(it->second).second; + full_array[index]=face_list[i]; + index++; + } + + } +} + +//!converts a list of global cell numbers +//!to a distributed array with local cell numbers +void ParallelTopology::convertGlobalCellList(const int* cell_list, int nbcell, int* local, int* ip) +{ + for (int i=0; i< nbcell; i++) + { + hash_map >::const_iterator iter = m_glob_to_loc.find(cell_list[i]); + ip[i]=(iter->second).first; + local[i]=(iter->second).second; + } +} + +/*!Converts a list of global face numbers + * to a distributed array with local face numbers + */ +void ParallelTopology::convertGlobalFaceList(const int* face_list, int nbface, int* local, int* ip) +{ + for (int i=0; i< nbface; i++) + { + hash_map >::const_iterator iter = m_face_glob_to_loc.find(face_list[i]); + if (iter == m_face_glob_to_loc.end()) + { + throw MED_EXCEPTION("convertGlobalFaceList - Face not found"); + } + ip[i]=(iter->second).first; + local[i]=(iter->second).second; + // cout << " in convertGlobalFAceList face global "< ("< >::iterator mmiter; + pair range=m_face_glob_to_loc.equal_range(face_list[i]); + for (mmiter it=range.first; it !=range.second; it++) + { + int ipfound=(it->second).first; + if (ipfound==ip) + local[i]=(it->second).second; + + } + } +} + + +////creating node mapping +// void ParallelTopology::createNodeMapping(map& type_connectivity, +// map& present_type_numbers, +// vector& polygon_conn, +// vector& polygon_conn_index, +// vector& polyhedron_conn, +// vector& polyhedron_conn_index, +// vector& polyhedron_face_index, +// int idomain) +// { +// set local_numbers; +// int local_index=0; + +// map::const_iterator iter; +// for (iter = present_type_numbers.begin(); iter!=present_type_numbers.end();iter++) +// { +// int type=iter->first; +// int nodes_per_type= type%100; + +// if (!((type/100==m_mesh_dimension) +// ||(type==MED_EN::MED_POLYGON && m_mesh_dimension==2) +// ||(type==MED_EN::MED_POLYHEDRA && m_mesh_dimension==3))) continue; + +// if (type != MED_EN::MED_POLYGON && type != MED_EN::MED_POLYHEDRA) +// { +// for (int icell=0; icell& globNodes) +// { +// // convert global nodes to local in the given domain +// set nodes; +// set::const_iterator n = globNodes.begin(); +// for ( ; n != globNodes.end(); ++n ) +// nodes.insert( convertGlobalNode( *n, domain )); + +// const MED_EN::medConnectivity connType = MED_EN::MED_NODAL; +// const MED_EN::medEntityMesh entity = MED_EN::MED_CELL; + +// // loop on all types of cells +// const MEDMEM::MESH* mesh = new_collection.getMesh( domain ); +// int nbTypes = mesh->getNumberOfTypesWithPoly( entity ); +// MED_EN::medGeometryElement * types = mesh->getTypesWithPoly( entity ); +// for ( int t = 0; t < nbTypes; ++t ) +// { +// // get connectivity +// if ( !mesh->existConnectivityWithPoly( connType, entity )) +// continue; +// int nbCell = mesh->getNumberOfElementsWithPoly( entity, types[t] ); +// const int *conn, *index; +// switch ( types[t] ) +// { +// case MED_EN::MED_POLYGON: +// conn = mesh->getPolygonsConnectivity( connType, entity ); +// index = mesh->getPolygonsConnectivityIndex( connType, entity ); +// break; +// case MED_EN::MED_POLYHEDRA: +// conn = mesh->getPolyhedronConnectivity( connType ); +// index = mesh->getPolyhedronFacesIndex(); +// nbCell = mesh->getNumberOfPolyhedronFaces(); +// break; +// default: +// conn = mesh->getConnectivity(MED_EN::MED_FULL_INTERLACE,connType, entity, types[t]); +// index = mesh->getConnectivityIndex(connType, entity); +// } +// // find a cell containing the first of given nodes, +// // then check if the found cell contains all the given nodes +// const int firstNode = *nodes.begin(); +// for ( int i = 0; i < nbCell; ++i ) +// { +// for ( int j = index[i]-1; j < index[i+1]-1; ++j ) +// if ( conn[j] == firstNode ) +// { +// int nbSame = 0; +// for ( j = index[i]-1; j < index[i+1]-1; ++j ) +// nbSame += nodes.count( conn[j] ); +// if ( nbSame == nodes.size() ) +// return true; +// break; +// } +// } +// } +// delete [] types; +// return false; +// } + +////creating face mapping +// void ParallelTopology::createFaceMapping(const MESHCollection& initial_collection, +// const MESHCollection& new_collection) +// // map& type_list, +// // map& present_type_numbers, +// // int idomain + +// { +// // containers for the new topology +// vector new_counts(m_nb_domain,0); +// vector domain_counts(m_nb_domain,0); +// const Topology* old_topology=initial_collection.getTopology(); +// int nb_domain_old=old_topology->nbDomain(); +// int global_index=old_topology->getFaceNumber(); +// //cout << "nb faces = " << global_index << endl; +// set > > global_treated; + +// //definition of the d-1 constituent for the considered mesh dimension +// MED_EN::medEntityMesh constituent_entity; +// switch (m_mesh_dimension) +// { +// case 3: +// constituent_entity= MED_EN::MED_FACE; +// break; +// case 2: +// constituent_entity = MED_EN::MED_EDGE; +// break; +// } + +// for (int iold=0; ioldgetNumberOfElements(constituent_entity,MED_EN::MED_ALL_ELEMENTS); +// int nbtotalface = initial_collection.getMesh(iold)->getNumberOfElementsWithPoly(constituent_entity,MED_EN::MED_ALL_ELEMENTS); +// SCRUTE_MED(nbplainface); +// SCRUTE_MED(nbtotalface); +// const int* face_conn; +// const int* face_offset; +// const int* poly_conn; +// const int* poly_index; +// if (nbtotalface >0) +// { +// if (nbplainface >0) +// { +// face_conn = initial_collection.getMesh(iold)->getConnectivity(MED_EN::MED_FULL_INTERLACE, +// MED_EN::MED_NODAL,constituent_entity,MED_EN::MED_ALL_ELEMENTS); +// face_offset = initial_collection.getMesh(iold)->getConnectivityIndex(MED_EN::MED_NODAL,constituent_entity); +// } +// if (nbtotalface > nbplainface) +// { +// poly_conn = initial_collection.getMesh(iold)->getPolygonsConnectivity(MED_EN::MED_NODAL,constituent_entity); +// poly_index = initial_collection.getMesh(iold)->getPolygonsConnectivityIndex(MED_EN::MED_NODAL,constituent_entity); +// } + +// } +// else +// { +// face_conn=0; +// face_offset=0; +// } +// for (int iface=0;ifaceconvertFaceToGlobal(iold,iface+1); + +// // int inode = face_offset[iface]; +// for (int i=0; i nodes; +// int nbnodes; +// if (ifaceconvertNodeToGlobal(iold,node); +// // cout << "global node "< >::iterator mmiter; +// pair range=m_node_glob_to_loc.equal_range(global); + +// int ip; +// for (mmiter it=range.first; it !=range.second; it++) +// { +// ip=(it->second).first; +// domain_counts[ip]++; +// } +// } +// } +// else +// { +// nbnodes =poly_index[iface-nbplainface+1]- poly_index[iface-nbplainface]; +// for (int inode= poly_index[iface-nbplainface];inode < poly_index[iface-nbplainface+1]; inode++) +// { +// int node=poly_conn[inode-1]; +// // SCRUTE_MED(node); +// int global = old_topology->convertNodeToGlobal(iold,node); +// // cout << "global node "< >::iterator mmiter; +// pair range=m_node_glob_to_loc.equal_range(global); + +// int ip; +// for (mmiter it=range.first; it !=range.second; it++) +// { +// ip=(it->second).first; +// domain_counts[ip]++; +// } +// } +// } +// set::const_iterator iter_node = nodes.begin(); +// int numbers[3]; +// numbers[2]=0; // for segments +// for (int i=0; i > > ::iterator iter_triplets; +// pair > triplet = make_pair(numbers[0],make_pair(numbers[1],numbers[2])); +// iter_triplets=global_treated.find(triplet); +// if (iter_triplets==global_treated.end()) +// { +// global_treated.insert(triplet); +// // int nbnodes=face_offset[iface+1]-face_offset[iface]; +// if (global_face_number == -1) +// { +// global_index++; +// global_face_number=global_index; + +// } +// // SCRUTE_MED(nbnodes); + +// for (int inew=0;inew new_counts(m_nb_domain,0); +// vector domain_counts(m_nb_domain,0); +// const Topology* old_topology=initial_collection.getTopology(); +// int nb_domain_old=old_topology->nbDomain(); +// //int global_index=old_topology->getFaceNumber(); +// // set > > global_treated; + +// //definition of the d-1 constituent for the considered mesh dimension +// MED_EN::medEntityMesh constituent_entity; +// switch (m_mesh_dimension) +// { +// case 3: +// constituent_entity= MED_EN::MED_FACE; +// break; +// case 2: +// constituent_entity = MED_EN::MED_EDGE; +// break; +// } + +// for (int iold=0; ioldgetCellNumber(iold); + +// const int* face_conn = initial_collection.getMesh(iold)-> +// getConnectivity(MED_EN::MED_FULL_INTERLACE, +// MED_EN::MED_DESCENDING,MED_EN::MED_CELL,MED_EN::MED_ALL_ELEMENTS); +// const int* face_offset = initial_collection.getMesh(iold)-> +// getConnectivityIndex(MED_EN::MED_DESCENDING,MED_EN::MED_CELL); +// MESSAGE_MED("end of connectivity calculation"); +// set global_treated; +// for (int icell=0; icellconvertCellToGlobal(iold,icell+1); +// int inew = getCellDomainNumber(global_cell_number); + +// for (int iface = face_offset[icell]; iface < face_offset[icell+1]; iface++) +// { +// int global_face_number=old_topology->convertFaceToGlobal(iold,abs(face_conn[iface-1])); +// if (global_treated.find (global_face_number)==global_treated.end()) +// { +// new_counts[inew]++; +// m_face_glob_to_loc.insert(make_pair(global_face_number,make_pair(inew,new_counts[inew]))); + +// //m_face_loc_to_glob.insert(make_pair(make_pair(inew,new_counts[inew]),global_face_number)); +// m_face_loc_to_glob[inew].push_back(global_face_number); +// global_treated.insert(global_face_number); + +// } +// } +// } +// } + + +// for (int inew=0;inew& type_connectivity, +// map& present_type_numbers, +// int idomain, +// MED_EN::medEntityMesh entity) +// { +// int dimension; +// switch (entity) +// { +// case MED_EN::MED_CELL: +// dimension=m_mesh_dimension; +// break; +// case MED_EN::MED_FACE: +// dimension=2; +// break; +// case MED_EN::MED_EDGE: +// dimension=1; +// break; +// } + +// MED_EN::MESH_ENTITIES::const_iterator currentEntity; +// list::const_iterator iter; +// currentEntity = MED_EN::meshEntities.find(MED_EN::MED_CELL); + +// for (iter = (*currentEntity).second.begin();iter != (*currentEntity).second.end(); iter++) +// { +// MED_EN::medGeometryElement type = (*iter); +// if (type/100 != dimension) continue; +// for (int inode=0; inodegetGlobal()<<" - "<getGlobal(),make_pair(idomain,ilocal))); +// m_face_loc_to_glob[idomain].push_back(face->getGlobal()); +// ilocal++; +// } +// } +// m_nb_faces[idomain] =ilocal-1; +// } +// } + +//================================================================================ +/*! + * \brief Recreating cell and node mapping after send-reveive and fusion of domain meshes + */ +//================================================================================ + +// void ParallelTopology::recreateMappingAfterFusion(const vector& meshes) +// { +// const char* LOC = "ParallelTopology::recreateMappingAfterFusion(): "; + +// m_glob_to_loc.clear(); +// m_node_glob_to_loc.clear(); +// m_face_glob_to_loc.clear(); + +// for (int idomain=0; idomaingetCoordinateptr() ) continue; // empty domian + +// //creating cell maps + +// m_nb_cells[idomain]=meshes[idomain]->getNumberOfElementsWithPoly(MED_EN::MED_CELL, +// MED_EN::MED_ALL_ELEMENTS); +// if ( m_cell_loc_to_glob_fuse[idomain].size() != m_nb_cells[idomain] ) +// throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused cells"); + +// m_loc_to_glob[idomain].swap( m_cell_loc_to_glob_fuse[idomain] ); + +// for (int i=0; igetNumberOfNodes(); +// m_node_loc_to_glob[idomain] = ((MEDMEM::MeshFuse*)meshes[idomain])->getNodeNumbers(); +// if ( m_node_loc_to_glob[idomain].size() != m_nb_nodes[idomain] ) +// throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused nodes"); + +// // setting mappings for all nodes +// for (int inode=0; inodegetNumberOfElementsWithPoly(constituent_entity, +// MED_EN::MED_ALL_ELEMENTS); +// if ( m_face_loc_to_glob_fuse[idomain].size() != m_nb_faces[idomain] ) +// throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused faces of domain "<< idomain +// << ": expect " << m_nb_faces[idomain] +// << " but have " << m_face_loc_to_glob_fuse[idomain].size()); + +// m_face_loc_to_glob[idomain].swap( m_face_loc_to_glob_fuse[idomain] ); + +// for (int iface=0; iface +#include +#include +//#include "boost/shared_ptr.hpp" + +#include "MEDPARTITIONER_Topology.hxx" + +namespace __gnu_cxx +{ + template<> struct hash< pair > + { + size_t operator()( const pair& x ) const + { + return hash< int >()( x.first*1000000+x.second ); + } + }; +} + +using namespace std; +namespace MEDPARTITIONER { + + class Graph; + class MESHCollection; + class MEDPARTITIONER_FaceModel; + + class ParallelTopology:public Topology + { + + public: + + ParallelTopology(); + + ParallelTopology(const std::vector&, + const std::vector&, + std::vector&, + std::vector&, + std::vector&); + + ParallelTopology(Graph* graph, int nbdomain, int mesh_dimension); + + ~ParallelTopology(); + //!converts a list of global cell numbers + //!to a distributed array with local cell numbers + void convertGlobalNodeList(const int*, int,int*,int*); + void convertGlobalNodeList(const int*, int,int*,int); + void convertGlobalNodeListWithTwins(const int* face_list, int nbnode, int*& local, int*& ip, int*& full_array, int& size); + + //!converts a list of global node numbers + //!to a distributed array with local cell numbers + void convertGlobalCellList(const int*, int , int*, int *); + + //!converts a list of global face numbers + //!to a distributed array with local face numbers + void convertGlobalFaceList(const int*, int , int*, int *); + void convertGlobalFaceList(const int*, int , int*, int); + void convertGlobalFaceListWithTwins(const int* face_list, int nbface, int*& local, int*& ip, int*& full_array,int& size); + + //!creating node mapping + +// void createNodeMapping(std::map& type_connectivity, +// std::map& present_type_numbers, +// std::vector& conn, +// std::vector& conn_index, +// std::vector& polyhedron_conn, +// std::vector& polyhedron_conn_index, +// std::vector& polyhedron_face_index, +// int idomain); + + // void createFaceMapping(std::map& type_connectivity, + // std::map& present_type_numbers, + // int idomain); + +// void createFaceMapping(const MESHCollection &, const MESHCollection&); +// void createFaceMapping2ndversion(const MESHCollection &); + +// //!converting node global numberings to local numberings +// void convertToLocal(std::map& type_connectivity, +// std::map& present_type_numbers, +// int idomain, +// MED_EN::medEntityMesh entity); + void convertToLocal2ndVersion(int* nodes, int nbnodes, int idomain); + + + + //! computing arrays with node/node correspondencies + // void computeNodeNodeCorrespondencies(int nbdomain,std::vector& ) const; + + //! computing arrays with node/node correspondencies + // void computeCellCellCorrespondencies(int nbdomain,std::vector&, const Graph*) const; + + //! retrieving Graph + // boost::shared_ptr getGraph() const; + + + //!converting node local numbering to global + inline int convertNodeToGlobal(int ip,int icell) const + { + //return m_node_loc_to_glob.find(make_pair(ip,icell))->second; + return m_node_loc_to_glob[ip][icell]; + } + + //!converting face local numbering to global + inline int convertFaceToGlobal(int ip,int iface) const + { + // if (m_face_loc_to_glob.find(make_pair(ip,icell))==m_face_loc_to_glob.end()) + // return -1; + // else + //return m_face_loc_to_glob.find(make_pair(ip,icell))->second; + return m_face_loc_to_glob[ip][iface]; + } + + //converting cell global numbering to local + inline int convertCellToGlobal(int ip,int icell) const + { + // if (m_loc_to_glob.find(make_pair(ip,icell))==m_loc_to_glob.end()) + // return -1; + // else + //return m_loc_to_glob.find(make_pair(ip,icell))->second; + return m_loc_to_glob[ip][icell]; + } + + inline void convertNodeToGlobal(int ip, const int* local, int n, int* global)const + { + for (int i=0; isecond; + } + + inline void convertCellToGlobal(int ip, const int* local, int n, int* global)const + { + for (int i=0; isecond; + } + + inline void convertFaceToGlobal(int ip, const int* local, int n, int* global)const + { + for (int i=0; isecond; + } + } + + inline int nbDomain() const + { + return m_nb_domain; + } + + int nbCells() const + { + return m_nb_total_cells; + } + int nbNodes() const + {return m_nb_total_nodes;} + + inline int nbCells( int idomain) const + { + return m_nb_cells[idomain]; + } + + + + + //!retrieving number of nodes + inline int getNodeNumber(int idomain) const + { + return m_nb_nodes[idomain]; + } + + inline int getNodeNumber() const + { + if (m_node_glob_to_loc.empty()) return 0; + set keys; + for (__gnu_cxx::hash_multimap >::const_iterator iter= m_node_glob_to_loc.begin(); + iter!=m_node_glob_to_loc.end(); + iter++) { + keys.insert(iter->first); + } + return keys.size(); + } + + //!retrieving list of nodes in global numbers + inline void getNodeList(int idomain, int* list) const + { + for (int i=0; isecond; + } + } + + //!< retrieving cell numbers after fusing in parallel mode + std::vector & getFusedCellNumbers(int idomain) + { + return m_cell_loc_to_glob_fuse[idomain]; + } + const std::vector & getFusedCellNumbers(int idomain) const + { + return m_cell_loc_to_glob_fuse[idomain]; + } + + //!< retrieving face numbers after fusing in parallel mode + std::vector & getFusedFaceNumbers(int idomain) + { + return m_face_loc_to_glob_fuse[idomain]; + } + const std::vector & getFusedFaceNumbers(int idomain) const + { + return m_face_loc_to_glob_fuse[idomain]; + } + + + //!retrieving number of nodes + inline int getCellNumber(int idomain) const + { + return m_nb_cells[idomain]; + } + + inline int getCellDomainNumber(int global) const + { + return (m_glob_to_loc.find(global)->second).first; + } + + //!retrieving list of nodes in global numbers + inline void getCellList(int idomain, int* list) const + { + for (int i=0; isecond; + } + + } + + inline int getFaceNumber(int idomain) const + { + return m_nb_faces[idomain]; + } + + inline int getFaceNumber() const + { + if (m_face_glob_to_loc.empty()) return 0; + set keys; + for (__gnu_cxx::hash_multimap >::const_iterator iter= m_face_glob_to_loc.begin(); + iter!=m_face_glob_to_loc.end(); + iter++) { + keys.insert(iter->first); + } + return keys.size(); + } + + + //!retrieving list of faces in global numbers + inline void getFaceList(int idomain, int* list) const + { + for (int i=0; isecond; + } + + } + + //! converting a global cell number to a local representation (domain + local number) + inline std::pair convertGlobalCell(int iglobal) const + { + return m_glob_to_loc.find(iglobal)->second; + } + + inline int convertGlobalFace(int iglobal, int idomain) + { + typedef __gnu_cxx::hash_multimap >::const_iterator MMiter; + std::pair eq = m_face_glob_to_loc.equal_range(iglobal); + for (MMiter it=eq.first; it != eq.second; it++) + { + // SCRUTE_MED (it->second.first); + //SCRUTE_MED (idomain); + if (it->second.first == idomain) return it->second.second; + + } + return -1; + } + + inline int convertGlobalNode(int iglobal, int idomain) + { + typedef __gnu_cxx::hash_multimap >::const_iterator MMiter; + pair eq = m_node_glob_to_loc.equal_range(iglobal); + for (MMiter it=eq.first; it != eq.second; it++) + { + if (it->second.first == idomain) return it->second.second; + } + return -1; + } + //!adding a face to the topology + inline void appendFace(int idomain, int ilocal, int iglobal) + { + m_face_loc_to_glob[idomain].push_back(iglobal); + m_face_glob_to_loc.insert(make_pair(iglobal,make_pair(idomain,ilocal))); + } + + //return max global face number + int getMaxGlobalFace() const; + + + //!recreating a face mapping from scratch + // void recreateFaceMapping(const TGeom2FacesByDomain& face_map); + + // recreating cell and node mapping after send-reveive and fusion of domain meshes + // virtual void recreateMappingAfterFusion(const std::vector& ); + + + + private: + + bool hasCellWithNodes( const MESHCollection&, int dom, const std::set& nodes ); + + + private: + //!mapping global -> local + typedef __gnu_cxx::hash_map > TGlob2DomainLoc; + TGlob2DomainLoc m_glob_to_loc; + + // bool is_equal_pair (pair a, pair b){ + // return (a.first==b.first && a.second==b.second); + // } + //!mapping local -> global + //map,int> m_loc_to_glob; + + // + //__gnu_cxx::hash_map,int, __gnu_cxx::hash > > m_loc_to_glob; + std::vector > m_loc_to_glob; + //!mapping global -> local + __gnu_cxx::hash_multimap > m_node_glob_to_loc; + + //!mapping local -> global + // map,int> m_node_loc_to_glob; + //__gnu_cxx::hash_map,int, __gnu_cxx::hash > > m_node_loc_to_glob; + std::vector > m_node_loc_to_glob; + + // global numbers in parallel mode + std::vector > m_cell_loc_to_glob_fuse; // glob nums after fusing + std::vector > m_face_loc_to_glob_fuse; // glob nums after fusing + + + //!mapping global -> local + typedef __gnu_cxx::hash_multimap > TGlob2LocsMap; + TGlob2LocsMap m_face_glob_to_loc; + + //!mapping local -> global + //__gnu_cxx::hash_map,int, __gnu_cxx::hash > > m_face_loc_to_glob; + std::vector > m_face_loc_to_glob; + + //map,int> m_face_loc_to_glob; + + std::vector m_nb_cells; + + std::vector m_nb_nodes; + + std::vector m_nb_faces; + + int m_nb_total_cells; + + int m_nb_total_nodes; + + int m_nb_total_faces; + + int m_nb_domain; + + int m_mesh_dimension; + + }; + + +} +#endif /*PARALLELTOPOLOGY_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_SCOTCHGraph.cxx b/src/MEDPartitioner/MEDPARTITIONER_SCOTCHGraph.cxx new file mode 100644 index 000000000..9c6abb79a --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_SCOTCHGraph.cxx @@ -0,0 +1,105 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#include +extern "C" { +#define restrict +#include "bin/scotch.h" +} +#include "MEDPARTITIONER_Graph.hxx" +#include "MEDPARTITIONER_SCOTCHGraph.hxx" + +using namespace MEDPARTITIONER; + +SCOTCHGraph::SCOTCHGraph():Graph() +{ +} + +SCOTCHGraph::SCOTCHGraph(MEDPARTITIONER::MEDSKYLINEARRAY* graph, int* edgeweight):Graph(graph,edgeweight) +{ +} + +SCOTCHGraph::~SCOTCHGraph() +{ +} + +void SCOTCHGraph::partGraph(int ndomain, const std::string& options_string, ParaDomainSelector* sel) +{ + // number of graph vertices + int n = m_graph->getNumberOf(); + + //graph + int * xadj=const_cast(m_graph->getIndex()); + int * adjncy = const_cast(m_graph->getValue()); + + //ndomain + int nparts = ndomain; + + // output parameters + int* partition = new int[n+1]; + + SCOTCH_Graph scotch_graph; + + SCOTCH_graphInit(&scotch_graph); + + + SCOTCH_graphBuild(&scotch_graph, + 1, //premier indice 1 + n, // nb of graph nodes + xadj, + 0, + m_cellweight, //graph vertices loads + 0, + xadj[n], // number of edges + adjncy, + m_edgeweight); + + SCOTCH_Strat scotch_strategy; + SCOTCH_stratInit(&scotch_strategy); + + //!user-defined options for the strategy + if (options_string!="") + SCOTCH_stratGraphMap(&scotch_strategy,options_string.c_str()); + + + if (nparts>1) + SCOTCH_graphPart(&scotch_graph,nparts,&scotch_strategy,partition); + else + // partition for 1 subdomain + for (int i=0; i index(n+1); + std::vector value(n); + index[0]=0; + for (int i=0; i + +namespace MEDPARTITIONER { + class MEDSKYLINEARRAY; + class MEDPARTITIONER_EXPORT SCOTCHGraph:public Graph + { + public: + SCOTCHGraph(); + SCOTCHGraph(MEDPARTITIONER::MEDSKYLINEARRAY*, int* edgeweight=0); + virtual ~SCOTCHGraph(); + void partGraph(int ndomain, const std::string& options_string="", ParaDomainSelector* sel=0); + //private: + // const MEDMEM::MEDSKYLINEARRAY* m_graph; + }; +} +#endif /*SCOTCHGRAPH_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.cxx b/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.cxx new file mode 100644 index 000000000..3076ae4f5 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.cxx @@ -0,0 +1,82 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#include "MEDMEM_SkyLineArray.hxx" +#include "MEDMEM_ConnectZone.hxx" +#include "MEDMEM_DriversDef.hxx" +#include "MEDMEM_Mesh.hxx" +#include "MEDMEM_Exception.hxx" + +#include "MEDPARTITIONER_Topology.hxx" +#include "MEDPARTITIONER_SequentialTopology.hxx" + +using namespace MEDPARTITIONER; + +SequentialTopology::SequentialTopology(){throw MEDEXCEPTION("SequentialTopology - not implemented yet");} + +SequentialTopology::SequentialTopology(std::vector, std::vector){throw MEDEXCEPTION("SequentialTopology - not implemented yet");} + +//!converts a list of global cell numbers +//!to a distributed array with local cell numbers +void SequentialTopology::convertGlobalNodeList(const int*, int,int*,int*){} + +//!converts a list of global node numbers +//!to a distributed array with local cell numbers +void SequentialTopology::convertGlobalCellList(const int*, int , int*, int *){} + +//number of doamins +int SequentialTopology::nbDomain() const +{ + return 0; +} + +//number of cells +int SequentialTopology::nbCells() const +{ + return 0; +} + +//number of cells on a specific domain +int SequentialTopology::nbCells(int idomain) const +{ + return 0; +} + +//!creating node mapping +void SequentialTopology::createNodeMapping(vector type_connectivity,int* present_type_numbers, int idomain){} + +//!converting node global numberings to local numberings +void SequentialTopology::convertToLocal(vector type_connectivity,int* present_type_numbers){} + +//!retrieving number of nodes +int SequentialTopology::getNodeNumber(int idomain) const +{ + return 0; +} + +//!retrieving list of nodes +void SequentialTopology::getNodeList(int idomain, int* list) const{} + +//!retrieving number of cells +int SequentialTopology::getCellNumber(int idomain) const +{ + return 0; +} + +//!retrieving list of cells +void SequentialTopology::getCellList(int idomain, int* list) const{} diff --git a/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.hxx b/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.hxx new file mode 100644 index 000000000..765819b73 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_SequentialTopology.hxx @@ -0,0 +1,90 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef SequentialTOPOLOGY_HXX_ +#define SequentialTOPOLOGY_HXX_ + +#include "MEDPARTITIONER.hxx" + +namespace MEDPARTITIONER { + class SequentialTopology:public Topology + { + + public: + + SequentialTopology(); + SequentialTopology(std::vector, std::vector); + + //!converts a list of global cell numbers + //!to a distributed array with local cell numbers + void convertGlobalNodeList(const int*, int,int*,int*); + + //!converts a list of global node numbers + //!to a distributed array with local cell numbers + void convertGlobalCellList(const int*, int , int*, int *); + + //number of doamins + int nbDomain() const; + + //number of cells + int nbCells() const; + + //number of cells on a specific domain + int nbCells(int idomain) const; + + //!creating node mapping + void createNodeMapping(vector type_connectivity,int* present_type_numbers, int idomain); + + //!converting node global numberings to local numberings + void convertToLocal(vector type_connectivity,int* present_type_numbers); + + //!retrieving number of nodes + int getNodeNumber(int idomain) const ; + + //!retrieving list of nodes + void getNodeList(int idomain, int* list) const; + + //!retrieving number of cells + int getCellNumber(int idomain) const ; + + //!retrieving list of cells + void getCellList(int idomain, int* list) const; + + private: + //!mapping global -> local + map > m_glob_to_loc; + + //!mapping local -> global + map,int> m_loc_to_glob; + + //!mapping global -> local + multimap > m_node_glob_to_loc; + + //!mapping local -> global + map,int> m_node_loc_to_glob; + + vector m_nb_cells; + + vector m_nb_nodes; + + int m_nb_total_cells; + + int m_nb_total_nodes; + }; +} +#endif /*SequentialTOPOLOGY_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.cxx b/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.cxx new file mode 100644 index 000000000..20cfd4c50 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.cxx @@ -0,0 +1,56 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +#include "MEDPARTITIONER_SkyLineArray.hxx" +#include "MEDMEM_Utilities.hxx" +#include + + +using namespace MEDPARTITIONER; + +MEDSKYLINEARRAY::MEDSKYLINEARRAY() +{ + MESSAGE_MED("Constructeur MEDSKYLINEARRAY sans parametre"); +} + +MEDSKYLINEARRAY::MEDSKYLINEARRAY(const MEDSKYLINEARRAY &myArray) +{ + _index=myArray._index; + _value=myArray._value; +} + +MEDSKYLINEARRAY::~MEDSKYLINEARRAY() +{ + MESSAGE_MED("Destructeur ~MEDSKYLINEARRAY"); + + //if (_index != NULL) delete [] _index; + //if (_value != NULL) delete [] _value; +} + + +MEDSKYLINEARRAY::MEDSKYLINEARRAY(const std::vector& index, const std::vector& value) +{ + _value=value; + _index=index; +} + + diff --git a/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.hxx b/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.hxx new file mode 100644 index 000000000..21d6fffbc --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_SkyLineArray.hxx @@ -0,0 +1,83 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// + +# ifndef __MEDPARTITIONER_MEDSKYLINEARRAY_H__ +# define __MEDPARTITIONER_MEDSKYLINEARRAY_H__ + +#include "MEDPARTITIONER.hxx" + +#include + +namespace MEDPARTITIONER { +class MEDPARTITIONER_EXPORT MEDSKYLINEARRAY +{ +private : + std::vector _index; + std::vector _value; + +public : + // Attention, avec ce constructeur, il n'est possible de remplir le MEDSKYLINEARRAY + MEDSKYLINEARRAY(); + + // Constructeur par recopie + MEDSKYLINEARRAY( const MEDSKYLINEARRAY &myArray ); + + // Avec ce constructeur la mémoire pour le tableau de valeur et le + // tableau d'index est réservée. Il suffit d'effectuer les séquences + // d'appels suivantes pour initialiser le MEDSKYLINEARRAY + // 1) setIndex(index) puis fois setI(i,&listValeurN°I) avec i dans 1..count + // rem : listValeurN°I est dupliquée + // 2) appeler fois setIJ(i,j,valeur) avec i dans 1..count et avec j dans 1..count + MEDSKYLINEARRAY( const std::vector& index, const std::vector& value ); + + + + ~MEDSKYLINEARRAY(); + //void setMEDSKYLINEARRAY( const int count, const int length, int* index , int* value ) ; + + inline int getNumberOf() const; + inline int getLength() const; + inline const int* getIndex() const; + inline const int* getValue() const; +}; + +// --------------------------------------- +// Methodes Inline +// --------------------------------------- +inline int MEDSKYLINEARRAY::getNumberOf() const +{ + return _index.size()-1; +} +inline int MEDSKYLINEARRAY::getLength() const +{ + return _value.size() ; +} +inline const int* MEDSKYLINEARRAY::getIndex() const +{ + return (const int*)(&_index[0]) ; +} +inline const int* MEDSKYLINEARRAY::getValue() const +{ + return (const int*)(&_value[0]) ; +} +} +# endif diff --git a/src/MEDPartitioner/MEDPARTITIONER_Topology.hxx b/src/MEDPartitioner/MEDPARTITIONER_Topology.hxx new file mode 100644 index 000000000..a38bdda02 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_Topology.hxx @@ -0,0 +1,173 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MEDPARTITIONER_TOPOLOGY_HXX_ +#define MEDPARTITIONER_TOPOLOGY_HXX_ + +//#include "MEDMEM_define.hxx" + +//#include "boost/shared_ptr.hpp" + +#include +#include + +namespace MEDPARTITIONER +{ + class CONNECTZONE; + class MEDSKYLINEARRAY; +} +namespace ParaMEDMEM { + class MEDCouplingUMesh; +} + +namespace MEDPARTITIONER { + + class Graph; + class MESHCollection; + class MEDPARTITIONER_FaceModel; + + // typedef std::map > TGeom2Faces; + // typedef std::vector< TGeom2Faces > TGeom2FacesByDomain; + + class Topology + { + public: + Topology(){} + Topology(std::vector, std::vector){} + + virtual ~Topology(){} + + //!converts a list of global cell numbers + //!to a distributed array with local cell numbers + virtual void convertGlobalNodeList(const int* list, int nb, int* local, int*ip)=0; + virtual void convertGlobalNodeList(const int* list, int nb, int* local, int ip)=0; + //!converts a list of global node numbers + //!to a distributed array with local cell numbers + virtual void convertGlobalCellList(const int*list , int nb, int* local, int*ip)=0; + + //!converts a list of global face numbers + //!to a distributed array with local face numbers + virtual void convertGlobalFaceList(const int*list , int nb, int* local, int*ip)=0; + virtual void convertGlobalFaceList(const int*list , int nb, int* local, int ip)=0; + virtual void convertGlobalFaceListWithTwins(const int* face_list, int nbface, int*& local, int*& ip, int*& full_array, int& size)=0; + virtual void convertGlobalNodeListWithTwins(const int* face_list, int nbnode, int*& local, int*& ip, int*& full_array, int& size)=0; + + + //number of doamins + virtual int nbDomain() const =0; + + //number of cells + virtual int nbCells() const=0; + + //number of nodes + virtual int nbNodes() const=0; + + //number of cells on a specific domain + virtual int nbCells(int idomain) const=0; + + // ////creating node mapping + // virtual void createNodeMapping(std::map& type_connectivity, + // std::map& present_type_numbers, + // std::vector& polygon_conn, +// std::vector& polygon_conn_index, +// std::vector& polyhedron_conn, +// std::vector& polyhedron_conn_index, +// std::vector& polyhedron_face_index, +// int domain)=0; + + ////creating face mapping + // virtual void createFaceMapping(std::map& type_connectivity, + // std::map& present_type_numbers,int domain)=0; + // + // virtual void createFaceMapping(const MESHCollection&,const MESHCollection&)=0; + + // //converting node global numberings to local numberings + // virtual void convertToLocal(std::map& type_connectivity, + // std::map& present_type_numbers, + // int idomain, + // MED_EN::medEntityMesh entity)=0; + //converting node global numberings to local numberings + virtual void convertToLocal2ndVersion(int*,int,int)=0; + + virtual int convertNodeToGlobal(int ip,int icell)const=0; + virtual int convertFaceToGlobal(int ip,int icell)const=0; + virtual int convertCellToGlobal(int ip,int icell)const=0; + + virtual void convertNodeToGlobal(int ip,const int* local, int n, int* global)const=0 ; + virtual void convertCellToGlobal(int ip,const int* local, int n, int* global)const=0 ; + virtual void convertFaceToGlobal(int ip,const int* local, int n, int* global)const=0 ; + + //retrieving number of nodes + virtual int getNodeNumber(int idomain) const =0; + virtual int getNodeNumber() const=0; + //retrieving list of nodes + virtual void getNodeList(int idomain, int* list) const =0; + + virtual std::vector & getFusedCellNumbers(int idomain) = 0; + virtual const std::vector & getFusedCellNumbers(int idomain) const = 0; + + virtual std::vector & getFusedFaceNumbers(int idomain) = 0; + virtual const std::vector & getFusedFaceNumbers(int idomain) const = 0; + + //retrieving number of nodes + virtual int getCellNumber(int idomain) const =0; + + //retrieving list of nodes + virtual void getCellList(int idomain, int* list) const =0; + + //retrieving number of faces + virtual int getFaceNumber(int idomain) const =0; + virtual int getFaceNumber()const =0; + + //retrieving list of nodes + virtual void getFaceList(int idomain, int* list) const =0; + + //adding a face to the mapping + virtual void appendFace(int idomain, int ilocal, int iglobal)=0; + + //return max global face number + virtual int getMaxGlobalFace()const=0; + + //return next free global face number + //virtual int nextGlobalFace(int start_num) const=0; + + //!converting a global cell number to a local representation + virtual std::pair convertGlobalCell(int iglobal) const =0; + + //converting a global face number to a local representation + virtual int convertGlobalFace(int iglobal, int idomain)=0; + + //converting a global node number to a local representation + virtual int convertGlobalNode(int iglobal, int idomain)=0; + + //! computing arrays with node/node correspondencies + // virtual void computeNodeNodeCorrespondencies(int nbdomain, std::vector&) const =0; + + //! computing arrays with cell/cell correspondencies + // virtual void computeCellCellCorrespondencies(int nbdomain, std::vector&, const Graph*) const =0; + + + //!recreating a face mapping from scratch + // virtual void recreateFaceMapping(const TGeom2FacesByDomain& )=0; + + //!recreating cell and node mapping after send-reveive and fusion of domain meshes + // virtual void recreateMappingAfterFusion(const std::vector& ) = 0; + }; +} + +#endif /*TOPOLOGY_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_UserGraph.cxx b/src/MEDPartitioner/MEDPARTITIONER_UserGraph.cxx new file mode 100644 index 000000000..9cbaeadea --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_UserGraph.cxx @@ -0,0 +1,54 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#include "MEDPARTITIONER_Graph.hxx" +#include "MEDPARTITIONER_UserGraph.hxx" +#include +#include +using namespace MEDPARTITIONER; + +/*! constructor that allows to specify the desired partition + * \param partition as table giving the domain number for each cell + * (domain numbers range from 0 to ndomain-1 + * \param n number of cells in the mesh + */ +UserGraph::UserGraph(MEDPARTITIONER::MEDSKYLINEARRAY* array, const int* partition, int n):Graph(array,0) +{ + + std::vector index(n+1),value(n); + + index[0]=0; + for (int i=0; i +namespace MEDPARTITIONER +{ + class MEDSKYLINEARRAY; + class ParaDomainSelector; + class MEDPARTITIONER_EXPORT UserGraph : public Graph + { + public: + UserGraph(MEDPARTITIONER::MEDSKYLINEARRAY*, const int*, int); + virtual ~UserGraph(); + void partGraph(int, const std::string& options=std::string(""), ParaDomainSelector* sel=0); + }; +} + +#endif /*MEDPARTITIONER_USERGRAPH_HXX_*/ diff --git a/src/MEDPartitioner/MEDPARTITIONER_utils.hxx b/src/MEDPartitioner/MEDPARTITIONER_utils.hxx new file mode 100644 index 000000000..1bb7d3069 --- /dev/null +++ b/src/MEDPartitioner/MEDPARTITIONER_utils.hxx @@ -0,0 +1,32 @@ +// Copyright (C) 2007-2010 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +#ifndef MEDPARTITIONER_UTILS_HXX_ +#define MEDPARTITIONER_UTILS_HXX_ + +#include + +namespace MEDPARTITIONER { + std::string trim(std::string& s,const std::string& drop = " ") + { + std::string r=s.erase(s.find_last_not_of(drop)+1); + return r.erase(0,r.find_first_not_of(drop)); + } + +} +#endif /*MEDPARTITIONER_UTILS_HXX_*/ diff --git a/src/MEDPartitioner/Makefile.am b/src/MEDPartitioner/Makefile.am new file mode 100644 index 000000000..8f47d930f --- /dev/null +++ b/src/MEDPartitioner/Makefile.am @@ -0,0 +1,130 @@ + # Copyright (C) 2007-2010 CEA/DEN, EDF R&D +# +# This library is free software; you can redistribute it and/or +# modify it under the terms of the GNU Lesser General Public +# License as published by the Free Software Foundation; either +# version 2.1 of the License. +# +# This library is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +# Lesser General Public License for more details. +# +# You should have received a copy of the GNU Lesser General Public +# License along with this library; if not, write to the Free Software +# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +# +# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +# +# MED MEDMEM : MED files in memory +# +include $(top_srcdir)/adm_local/unix/make_common_starter.am + +# this directory must be recompiled before Test folder + +if CPPUNIT_IS_OK + SUBDIRS=. +endif + +lib_LTLIBRARIES= libmedsplitter.la + +salomeinclude_HEADERS= \ +MEDPARTITIONER_Topology.hxx \ +MEDPARTITIONER_MESHCollection.hxx \ +MEDPARTITIONER_MESHCollectionDriver.H \ +MEDPARTITIONER_MESHCollectionMedXMLDriver.H \ +MEDPARTITIONER_MESHCollectionMedAsciiDriver.H \ +MEDPARTITIONER_MESHCollectionDriver.hxx \ +MEDPARTITIONER_MESHCollectionMedXMLDriver.hxx \ +MEDPARTITIONER_MESHCollectionMedAsciiDriver.hxx \ +MEDPARTITIONER_ParallelTopology.hxx \ +MEDPARTITIONER_FaceModel.hxx \ +MEDPARTITIONER_Graph.hxx\ +MEDPARTITIONER_UserGraph.hxx\ +MEDPARTITIONER_SequentialTopology.hxx \ +MEDPARTITIONER_utils.hxx \ +MEDPARTITIONER.hxx \ +MEDPARTITIONER_ParaDomainSelector.hxx \ +MEDPARTITIONER_ConnectZone.hxx \ +MEDPARTITIONER_SkyLineArray.hxx + +if MED_ENABLE_METIS + salomeinclude_HEADERS+= MEDPARTITIONER_METISGraph.hxx +endif +if MED_ENABLE_SCOTCH + salomeinclude_HEADERS+= MEDPARTITIONER_SCOTCHGraph.hxx +endif + +dist_libmedsplitter_la_SOURCES= \ +MEDPARTITIONER_MESHCollection.cxx \ +MEDPARTITIONER_MESHCollectionDriver.cxx \ +MEDPARTITIONER_MESHCollectionMedXMLDriver.cxx \ +MEDPARTITIONER_MESHCollectionMedAsciiDriver.cxx \ +MEDPARTITIONER_ParallelTopology.cxx \ +MEDPARTITIONER_Graph.cxx\ +MEDPARTITIONER_UserGraph.cxx\ +MEDPARTITIONER_ParaDomainSelector.cxx \ +MEDPARTITIONER_JointExchangeData.cxx \ +MEDPARTITIONER_SkyLineArray.cxx \ +MEDPARTITIONER_ConnectZone.cxx + +if MED_ENABLE_METIS + dist_libmedsplitter_la_SOURCES+= MEDPARTITIONER_METISGraph.cxx +endif +if MED_ENABLE_SCOTCH + dist_libmedsplitter_la_SOURCES+= MEDPARTITIONER_SCOTCHGraph.cxx +endif + +libmedsplitter_la_CPPFLAGS= $(MPI_INCLUDES) $(MED2_INCLUDES) $(HDF5_INCLUDES) @CXXTMPDPTHFLAGS@ \ + $(BOOST_CPPFLAGS) $(LIBXML_INCLUDES) \ + -I$(srcdir)/../MEDMEM -I$(srcdir)/../MEDWrapper/V2_1/Core \ + -I$(srcdir)/../INTERP_KERNEL/Bases -I$(srcdir)/../MEDCoupling \ + -I$(srcdir)/../MEDLoader + +libmedsplitter_la_LDFLAGS= +#libmedsplitter_la_LDFLAGS= $(MED2_LIBS) $(HDF5_LIBS) $(STDLIB) $(LIBXML_LIBS) \ +# ../MEDMEM/libmedmem.la ../MEDWrapper/V2_1/Core/libmed_V2_1.la + +if MED_ENABLE_PARMETIS + libmedsplitter_la_CPPFLAGS+= $(PARMETIS_CPPFLAGS) + libmedsplitter_la_LDFLAGS+= $(PARMETIS_LIBS) +endif +if MED_ENABLE_METIS + libmedsplitter_la_CPPFLAGS+= $(METIS_CPPFLAGS) + libmedsplitter_la_LDFLAGS+= $(METIS_LIBS) +endif +if MED_ENABLE_SCOTCH + libmedsplitter_la_CPPFLAGS+= $(SCOTCH_CPPFLAGS) + libmedsplitter_la_LDFLAGS+= $(SCOTCH_LIBS) +endif +if MED_ENABLE_KERNEL + libmedsplitter_la_CPPFLAGS+= ${KERNEL_CXXFLAGS} + libmedsplitter_la_LDFLAGS+= ${KERNEL_LDFLAGS} -lSALOMELocalTrace +endif + +libmedsplitter_la_LDFLAGS+= $(MED2_LIBS) $(HDF5_LIBS) $(STDLIB) $(LIBXML_LIBS) $(MPI_LIBS) \ + ../MEDMEM/libmedmem.la ../INTERP_KERNEL/libinterpkernel.la ../MEDCoupling/libmedcoupling.la ../MEDLoader/libmedloader.la + +# Executables targets +bin_PROGRAMS= medsplitter + +dist_medsplitter_SOURCES= medsplitter.cxx + +medsplitter_CPPFLAGS= $(libmedsplitter_la_CPPFLAGS) +medsplitter_LDADD= $(libmedsplitter_la_LDFLAGS) -lm $(BOOST_LIBS) libmedsplitter.la +if MED_ENABLE_KERNEL + medsplitter_LDADD+= -lSALOMEBasics +endif + +if MPI_IS_OK + bin_PROGRAMS+=medsplitter_para + dist_medsplitter_para_SOURCES= medsplitter_para.cxx + medsplitter_para_CPPFLAGS= $(medsplitter_CPPFLAGS) + medsplitter_para_LDADD= $(medsplitter_LDADD) +endif + +OBSOLETE_FILES = \ + MEDPARTITIONER_SequentialTopology.cxx \ + test_HighLevelAPI.cxx + +EXTRA_DIST += $(OBSOLETE_FILES)