Back to index

salome-med  6.5.0
MEDSPLITTER_ParallelTopology.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 #include <set>
00020 #include <map>
00021 #include <vector>
00022 
00023 #include "InterpKernelHashMap.hxx"
00024 
00025 #include "MEDMEM_CellModel.hxx"
00026 #include "MEDMEM_ConnectZone.hxx"
00027 #include "MEDMEM_DriversDef.hxx"
00028 #include "MEDMEM_Exception.hxx"
00029 #include "MEDMEM_Mesh.hxx"
00030 #include "MEDMEM_MeshFuse.hxx"
00031 #include "MEDMEM_SkyLineArray.hxx"
00032 #include "MEDMEM_Utilities.hxx"
00033 
00034 #include "MEDSPLITTER_MESHCollection.hxx"
00035 #include "MEDSPLITTER_Topology.hxx"
00036 #include "MEDSPLITTER_Graph.hxx"
00037 #include "MEDSPLITTER_ParallelTopology.hxx"
00038 
00039 using namespace INTERP_KERNEL;
00040 
00041 using namespace MEDSPLITTER;
00042 
00043 //empty constructor
00044 ParallelTopology::ParallelTopology():m_nb_domain(0),m_mesh_dimension(0)
00045 {}
00046 
00048 ParallelTopology::ParallelTopology(const vector<MEDMEM::MESH*>& meshes, 
00049                                    const vector<MEDMEM::CONNECTZONE*>& cz,
00050                                    vector<int*>& cellglobal,
00051                                    vector<int*>& nodeglobal,
00052                                    vector<int*>& faceglobal):m_nb_domain(meshes.size())/*,m_mesh_dimension(meshes[0]->getMeshDimension())*/
00053 {
00054 
00055   int index_global=0;
00056   int index_node_global=0;
00057   int index_face_global=0;
00058 
00059   m_nb_cells.resize(m_nb_domain);
00060   m_nb_nodes.resize(m_nb_domain);
00061   m_nb_faces.resize(m_nb_domain);
00062 
00063   m_loc_to_glob.resize(m_nb_domain);
00064   m_node_loc_to_glob.resize(m_nb_domain);
00065   m_face_loc_to_glob.resize(m_nb_domain);
00066 
00067   MED_EN::medEntityMesh constituent_entity;
00068 
00069   bool parallel_mode = false;
00070   for (int idomain=0; !parallel_mode && idomain<m_nb_domain; idomain++)
00071     parallel_mode = (!meshes[idomain]);
00072 
00073   for (int idomain=0; idomain<m_nb_domain; idomain++)
00074   {
00075     if ( !meshes[idomain] )
00076       continue;
00077     m_mesh_dimension = meshes[idomain]->getMeshDimension();
00078     constituent_entity = (m_mesh_dimension == 3 ? MED_EN::MED_FACE : MED_EN::MED_EDGE );
00079 
00080     //creating cell maps
00081     m_nb_cells[idomain]=meshes[idomain]->getNumberOfElements(MED_EN::MED_CELL, MED_EN::MED_ALL_ELEMENTS);
00082     //              cout << "Nb cells (domain "<<idomain<<") = "<<m_nb_cells[idomain];
00083     m_loc_to_glob[idomain].resize(m_nb_cells[idomain]);
00084 
00085     if (cellglobal[idomain]==0 || parallel_mode)
00086     {
00087       MESSAGE_MED("Creating global numbering"); 
00088       //creating global numbering from scratch
00089       for (int i=0; i<m_nb_cells[idomain]; i++)
00090       {
00091         index_global++;
00092         m_glob_to_loc[index_global]=make_pair(idomain,i+1);
00093         //m_loc_to_glob[make_pair(idomain,i+1)]=index_global;
00094         m_loc_to_glob[idomain][i]=index_global;
00095         //        cout<<"glob:"<<index_global<<" --> ("<<idomain<<","<<i+1<<")"<<endl;
00096       }
00097     }
00098     //using global numbering coming from a previous numbering
00099     else
00100     {
00101       MESSAGE_MED("Using former global numbering");
00102       for (int i=0; i<m_nb_cells[idomain]; i++)
00103       {
00104         int global=cellglobal[idomain][i];
00105         m_glob_to_loc[global]=make_pair(idomain,i+1);
00106         //m_loc_to_glob[make_pair(idomain,i+1)]=global;
00107         m_loc_to_glob[idomain][i]=global;
00108         index_global++;
00109         //        cout<<"glob:"<<global<<" --> ("<<idomain<<","<<i+1<<")"<<endl;
00110       }
00111     }
00112 
00113     //cas sequentiel
00114     if (m_nb_domain==1)
00115     {
00116       m_nb_total_cells=index_global;
00117       m_nb_cells[0]=index_global;
00118       m_node_loc_to_glob[idomain].resize(meshes[idomain]->getNumberOfNodes());
00119       for (int i=0; i<meshes[idomain]->getNumberOfNodes(); i++)
00120       {
00121         m_node_glob_to_loc.insert(make_pair(i+1,make_pair(0,i+1)));
00122         //m_node_loc_to_glob.insert(make_pair(make_pair(0,i+1), i+1));
00123         m_node_loc_to_glob[0][i]=i+1;
00124       }
00125       m_nb_total_nodes=meshes[idomain]->getNumberOfNodes();   
00126       m_nb_nodes[0]=m_nb_total_nodes; 
00127 
00128       //                      meshes[idomain]->getConnectivity( MED_EN::MED_DESCENDING, MED_EN::MED_CELL, MED_EN::MED_ALL_ELEMENTS); 
00129       int nbfaces=meshes[idomain]->getNumberOfElements(constituent_entity,MED_EN::MED_ALL_ELEMENTS);
00130       m_face_loc_to_glob[idomain].resize(nbfaces);
00131       for (int i=0; i<nbfaces; i++)
00132       {
00133         m_face_glob_to_loc.insert(make_pair(i+1,make_pair(0,i+1)));
00134         //m_face_loc_to_glob.insert(make_pair(make_pair(0,i+1), i+1));
00135         m_face_loc_to_glob[0][i]=i+1;
00136       }
00137       m_nb_total_faces=nbfaces;   
00138       m_nb_faces[0]=nbfaces;
00139       MESSAGE_MED ("nb total cells "<< m_nb_total_cells);
00140       MESSAGE_MED("nb total nodes "<< m_nb_total_nodes);  
00141       MESSAGE_MED("nb total faces "<< m_nb_total_faces);  
00142       return;
00143     }
00144 
00145     //creating node maps
00146     m_nb_nodes[idomain]=meshes[idomain]->getNumberOfNodes();
00147     INTERP_KERNEL::HashMap <int,pair<int,int> > local2distant;
00148     m_node_loc_to_glob[idomain].resize(m_nb_nodes[idomain]);
00149     for (unsigned icz=0; icz<cz.size(); icz++)
00150     {
00151       if (cz[icz]->getLocalDomainNumber() == idomain && 
00152           cz[icz]->getLocalDomainNumber()>cz[icz]->getDistantDomainNumber())
00153       {
00154         int nb_node= cz[icz]->getNodeNumber();
00155         const int* node_corresp=cz[icz]->getNodeCorrespValue();
00156         int distant_ip = cz[icz]->getDistantDomainNumber();
00157         for (int i=0; i< nb_node; i++)
00158         {
00159           int local= node_corresp[i*2];
00160           int distant = node_corresp[i*2+1];
00161           local2distant.insert(make_pair(local, make_pair(distant_ip,distant)));    
00162         }
00163       }
00164     }
00165     // setting mappings for all nodes
00166     if (nodeglobal[idomain]==0)
00167     {
00168       for (int inode=0; inode<m_nb_nodes[idomain]; inode++)
00169       {
00170         if (local2distant.find(inode+1)==local2distant.end())
00171         {
00172           index_node_global++;
00173           m_node_glob_to_loc.insert(make_pair(index_node_global,make_pair(idomain,inode+1)));
00174           //m_node_loc_to_glob[make_pair(idomain,inode+1)]=index_node_global;
00175           m_node_loc_to_glob[idomain][inode]=index_node_global;
00176         }   
00177         else
00178         {
00179           int ip = (local2distant.find(inode+1)->second).first;
00180           int distant = (local2distant.find(inode+1)->second).second;
00181           //int global_number=m_loc_to_glob[make_pair(ip,distant)];
00182           int global_number=m_loc_to_glob[ip][distant-1];
00183           m_node_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,inode+1)));
00184           //m_node_loc_to_glob[make_pair(idomain,inode+1)]=global_number;
00185           m_node_loc_to_glob[idomain][inode]=global_number;
00186         } 
00187       }
00188     }
00189     //using former node numbering
00190     else
00191     {//       cout << "("<<idomain<<","<<i+1<<")->"<<i+1<<endl;
00192       for (int inode=0; inode<m_nb_nodes[idomain]; inode++)
00193       {
00194         int global_number=nodeglobal[idomain][inode];
00195         //        cout << "global_number "<<global_number<<endl;
00196         m_node_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,inode+1)));
00197         //m_node_loc_to_glob[make_pair(idomain,inode+1)]=global_number;
00198         m_node_loc_to_glob[idomain][inode]=global_number;
00199       }
00200     }
00201 
00202 
00203     //creating  dimension d-1 component mappings
00204 
00205     //              meshes[idomain]->getConnectivity( MED_EN::MED_DESCENDING, MED_EN::MED_CELL, MED_EN::MED_ALL_ELEMENTS); 
00206     m_nb_faces[idomain]=meshes[idomain]->getNumberOfElements(constituent_entity,MED_EN::MED_ALL_ELEMENTS);
00207     MESSAGE_MED ("Nb faces domain " << idomain<<m_nb_faces[idomain]);
00208     m_face_loc_to_glob[idomain].resize(m_nb_faces[idomain]);
00209     local2distant.clear();
00210     for (unsigned icz=0; icz<cz.size(); icz++)
00211     {
00212       if (cz[icz]->getLocalDomainNumber() == idomain && 
00213           cz[icz]->getLocalDomainNumber()>cz[icz]->getDistantDomainNumber())
00214       {
00215         int nb_face= cz[icz]->getFaceNumber();
00216         const int* face_corresp=cz[icz]->getFaceCorrespValue();
00217         int distant_ip = cz[icz]->getDistantDomainNumber();
00218         for (int i=0; i< nb_face; i++)
00219         {
00220           int local= face_corresp[i*2];
00221           int distant = face_corresp[i*2+1];
00222           local2distant.insert(make_pair(local, make_pair(distant_ip,distant)));    
00223         }
00224       }
00225     }
00226     // setting mappings for all faces
00227     if (faceglobal[idomain]==0)
00228     {
00229       for (int iface=0; iface<m_nb_faces[idomain]; iface++)
00230       {
00231         if (local2distant.find(iface+1)==local2distant.end())
00232         {
00233           index_face_global++;
00234           m_face_glob_to_loc.insert(make_pair(index_face_global,make_pair(idomain,iface+1)));
00235           //m_face_loc_to_glob[make_pair(idomain,iface+1)]=index_face_global;
00236           m_face_loc_to_glob[idomain][iface]=index_face_global;
00237         }   
00238         else
00239         {
00240           int ip = (local2distant.find(iface+1)->second).first;
00241           int distant = (local2distant.find(iface+1)->second).second;
00242           //int global_number=m_loc_to_glob[make_pair(ip,distant)];
00243           int global_number=m_loc_to_glob[ip][distant-1];
00244           m_face_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,iface+1)));
00245           //m_face_loc_to_glob[make_pair(idomain,iface+1)]=global_number;
00246           m_face_loc_to_glob[idomain][iface]=global_number;
00247         } 
00248       }
00249     }
00250     //using former face numbering
00251     else
00252     {
00253       for (int iface=0; iface<m_nb_faces[idomain]; iface++)
00254       {
00255         int global_number=faceglobal[idomain][iface];
00256         //SCRUTE_MED(global_number);
00257         m_face_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,iface+1)));
00258         //m_face_loc_to_glob[make_pair(idomain,iface+1)]=global_number;
00259         m_face_loc_to_glob[idomain][iface]=global_number;
00260       }
00261     }
00262   }
00263 
00264   m_nb_total_cells=index_global;
00265   m_nb_total_nodes=index_node_global;   
00266   m_nb_total_faces=index_face_global;
00267   SCRUTE_MED(m_nb_total_cells);
00268   SCRUTE_MED(m_nb_total_faces);
00269   SCRUTE_MED(m_nb_total_nodes);
00270 
00271 }
00272 
00273 
00275 ParallelTopology::ParallelTopology(boost::shared_ptr<Graph> graph, int nb_domain, int mesh_dimension):
00276   m_nb_domain(nb_domain),
00277   m_mesh_dimension(mesh_dimension),
00278   m_nb_cells(graph->nbVertices()),
00279   m_graph(graph)
00280 {
00281   m_nb_cells.resize(m_nb_domain);
00282   m_nb_nodes.resize(m_nb_domain);
00283   m_nb_faces.resize(m_nb_domain);
00284 
00285   m_loc_to_glob.resize(m_nb_domain);
00286   m_node_loc_to_glob.resize(m_nb_domain);
00287   m_face_loc_to_glob.resize(m_nb_domain);
00288 
00289   // used in parallel mode only
00290   m_cell_loc_to_glob_fuse.resize(m_nb_domain);
00291   m_face_loc_to_glob_fuse.resize(m_nb_domain);
00292 
00293   for (int i=0; i<m_nb_domain; i++)
00294     m_nb_cells[i]=0;  
00295 
00296   const int* part = graph-> getPart();
00297   m_nb_total_cells= graph->nbVertices();
00298 
00299   for (int icell=0; icell<m_nb_total_cells; icell++)
00300   {
00301     int idomain = part[icell];
00302     m_nb_cells[idomain]++;
00303     //m_loc_to_glob[make_pair(idomain,m_nb_cells[idomain])]=icell+1;
00304     m_loc_to_glob[idomain].push_back(icell+1);
00305     m_glob_to_loc[icell+1]=make_pair(idomain,m_nb_cells[idomain]);
00306 
00307   }
00308   for (int idomain=0; idomain<m_nb_domain; idomain++)
00309     MESSAGE_MED("Nombre de cellules dans le domaine "<< idomain <<" : "<<m_nb_cells[idomain]); 
00310 
00311   SCRUTE_MED(m_nb_total_cells);
00312 
00313 }
00314 
00315 ParallelTopology::~ParallelTopology()
00316 {
00317 
00318 } 
00319 
00326 void ParallelTopology::convertGlobalNodeList(const int* node_list, int nbnode, int* local, int* ip)
00327 {
00328   if (m_node_glob_to_loc.empty()) 
00329     throw MEDMEM::MEDEXCEPTION("convertGlobalNodeList - Node mapping has not yet been built");
00330   for (int i=0; i< nbnode; i++)
00331   {
00332     pair<int,int> local_node = m_node_glob_to_loc.find(node_list[i])->second;
00333     ip[i]=local_node.first;
00334     local[i]=local_node.second;
00335   }
00336 }
00337 
00345 void ParallelTopology::convertGlobalNodeList(const int* node_list, int nbnode, int* local, int ip)
00346 {
00347   if (m_node_glob_to_loc.empty()) 
00348     throw MEDMEM::MEDEXCEPTION("convertGlobalNodeList - Node mapping has not yet been built");
00349 
00350   for (int i=0; i< nbnode; i++)
00351   {
00352     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00353     pair<mmiter,mmiter> range=m_node_glob_to_loc.equal_range(node_list[i]);
00354     for (mmiter it=range.first; it !=range.second; it++)
00355     { 
00356       int ipfound=(it->second).first;
00357       if (ipfound==ip)
00358         local[i]=(it->second).second;
00359     }
00360   }
00361 } 
00362 
00369 void ParallelTopology::convertGlobalNodeListWithTwins(const int* node_list, int nbnode, int*& local, int*& ip,int*& full_array, int& size)
00370 {
00371   if (m_node_glob_to_loc.empty()) 
00372     throw MEDMEM::MEDEXCEPTION("convertGlobalNodeList - Node mapping has not yet been built");
00373 
00374   size=0;
00375   for (int i=0; i< nbnode; i++)
00376   {
00377     int count= m_node_glob_to_loc.count(node_list[i]);
00378     //      if (count > 1) 
00379     //        cout << "noeud " << node_list[i]<< " doublon d'ordre " << count<<endl;
00380     size+=count;
00381   }
00382   int index=0;
00383   ip=new int[size];
00384   local=new int[size];
00385   full_array=new int[size];
00386   for (int i=0; i< nbnode; i++)
00387   {
00388     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00389     pair<mmiter,mmiter> range=m_node_glob_to_loc.equal_range(node_list[i]);
00390     for (mmiter it=range.first; it !=range.second; it++)
00391     { 
00392       ip[index]=(it->second).first;
00393       local[index]=(it->second).second;
00394       full_array [index]=node_list[i];
00395       index++;
00396     }
00397 
00398   }
00399 }
00400 
00407 void ParallelTopology::convertGlobalFaceListWithTwins(const int* face_list, int nbface, int*& local, int*& ip, int*& full_array,int& size)
00408 {
00409   size=0;
00410   for (int i=0; i< nbface; i++)
00411   {
00412     //int count = m_face_glob_to_loc.count(face_list[i]);
00413     //if (count >1) MESSAGE_MED("face en doublon "<<face_list[i]);
00414     size+= m_face_glob_to_loc.count(face_list[i]);
00415   }
00416   int index=0;
00417   ip=new int[size];
00418   local=new int[size];
00419   full_array=new int[size];
00420   for (int i=0; i< nbface; i++)
00421   {
00422     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00423     pair<mmiter,mmiter> range=m_face_glob_to_loc.equal_range(face_list[i]);
00424     for (mmiter it=range.first; it !=range.second; it++)
00425     { 
00426       ip[index]=(it->second).first;
00427       local[index]=(it->second).second;
00428       full_array[index]=face_list[i];
00429       index++;
00430     }
00431 
00432   }
00433 }
00434 
00437 void ParallelTopology::convertGlobalCellList(const int* cell_list, int nbcell, int* local, int* ip)
00438 {
00439   for (int i=0; i< nbcell; i++)
00440   {
00441     INTERP_KERNEL::HashMap<int, pair<int,int> >::const_iterator iter = m_glob_to_loc.find(cell_list[i]);
00442     ip[i]=(iter->second).first;
00443     local[i]=(iter->second).second;
00444   }
00445 }
00446 
00450 void ParallelTopology::convertGlobalFaceList(const int* face_list, int nbface, int* local, int* ip)
00451 {
00452   for (int i=0; i< nbface; i++)
00453   {
00454     INTERP_KERNEL::HashMap<int, pair<int,int> >::const_iterator iter = m_face_glob_to_loc.find(face_list[i]);
00455     if (iter == m_face_glob_to_loc.end())
00456     {
00457       throw MED_EXCEPTION("convertGlobalFaceList - Face  not found");
00458     }
00459     ip[i]=(iter->second).first;
00460     local[i]=(iter->second).second;
00461     //    cout << " in convertGlobalFAceList face global "<<face_list[i]<<" -> ("<<ip[i]<<","<<local[i]<<")"<<endl;
00462   }
00463 }
00464 
00472 void ParallelTopology::convertGlobalFaceList(const int* face_list, int nbface, int* local, int ip)
00473 {
00474   for (int i=0; i< nbface; i++)
00475   {
00476     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00477     pair<mmiter,mmiter> range=m_face_glob_to_loc.equal_range(face_list[i]);
00478     for (mmiter it=range.first; it !=range.second; it++)
00479     { 
00480       int ipfound=(it->second).first;
00481       if (ipfound==ip)
00482         local[i]=(it->second).second; 
00483 
00484     }
00485   }
00486 } 
00487 
00488 
00490 void ParallelTopology::createNodeMapping(map<MED_EN::medGeometryElement,int*>& type_connectivity,
00491                                          map<MED_EN::medGeometryElement,int>& present_type_numbers,
00492                                          vector<int>& polygon_conn,
00493                                          vector<int>& polygon_conn_index,
00494                                          vector<int>& polyhedron_conn,
00495                                          vector<int>& polyhedron_conn_index,
00496                                          vector<int>& polyhedron_face_index,
00497                                          int idomain)
00498 {
00499   set<int> local_numbers;
00500   int local_index=0;
00501 
00502   map<MED_EN::medGeometryElement,int>::const_iterator iter;
00503   for (iter = present_type_numbers.begin(); iter!=present_type_numbers.end();iter++)
00504   {
00505     int type=iter->first;
00506     int nodes_per_type= type%100;
00507 
00508     if (!((type/100==m_mesh_dimension)
00509           ||(type==MED_EN::MED_POLYGON && m_mesh_dimension==2)
00510           ||(type==MED_EN::MED_POLYHEDRA && m_mesh_dimension==3))) continue;
00511 
00512     if (type != MED_EN::MED_POLYGON && type != MED_EN::MED_POLYHEDRA)
00513     {
00514       for (int icell=0; icell<present_type_numbers[type]; icell++)
00515       {
00516         for (int inode=0; inode<nodes_per_type; inode++)
00517         {
00518           int global=type_connectivity[type][icell*nodes_per_type+inode];
00519           if(local_numbers.find(global)==local_numbers.end())
00520           {
00521             local_index++;
00522             local_numbers.insert(global);
00523             //m_node_loc_to_glob.insert(make_pair(make_pair(idomain,local_index),global));
00524             m_node_loc_to_glob[idomain].push_back(global);
00525             m_node_glob_to_loc.insert(make_pair(global,make_pair(idomain,local_index)));
00526             //          cout << "node : global ="<<global<<" local =("<<idomain<<","<<local_index<<")"<<endl;         
00527           }
00528         }
00529 
00530       }
00531     }
00532     else if (type== MED_EN::MED_POLYGON)
00533     {
00534       for ( unsigned i = 0; i < polygon_conn.size(); ++i )
00535       {
00536         int global=polygon_conn[i];
00537         if(local_numbers.find(global)==local_numbers.end())
00538         {
00539           local_index++;
00540           local_numbers.insert(global);
00541           m_node_loc_to_glob[idomain].push_back(global);
00542           m_node_glob_to_loc.insert(make_pair(global,make_pair(idomain,local_index)));
00543         }
00544       }
00545     }
00546     else if (type==MED_EN::MED_POLYHEDRA)
00547     {
00548       for ( unsigned i = 0; i < polyhedron_conn.size(); ++i )
00549       {
00550         int global=polyhedron_conn[i];
00551         if(local_numbers.find(global)==local_numbers.end())
00552         {
00553           local_index++;
00554           local_numbers.insert(global);
00555           m_node_loc_to_glob[idomain].push_back(global);
00556           m_node_glob_to_loc.insert(make_pair(global,make_pair(idomain,local_index)));
00557         }
00558       }
00559     }
00560 
00561   }
00562   m_nb_nodes[idomain]=local_index;
00563 }
00564 
00565 //================================================================================
00569 //================================================================================
00570 
00571 bool ParallelTopology::hasCellWithNodes( const MESHCollection& new_collection,
00572                                          int                   domain,
00573                                          const set<int>&       globNodes)
00574 {
00575   // convert global nodes to local in the given domain
00576   set<int> nodes;
00577   set<int>::const_iterator n = globNodes.begin();
00578   for ( ; n != globNodes.end(); ++n )
00579     nodes.insert( convertGlobalNode( *n, domain ));
00580 
00581   const MED_EN::medConnectivity connType = MED_EN::MED_NODAL;
00582   const MED_EN::medEntityMesh   entity   = MED_EN::MED_CELL;
00583 
00584   // loop on all types of cells
00585   const MEDMEM::MESH* mesh = new_collection.getMesh( domain );
00586   int nbTypes = mesh->getNumberOfTypes( entity );
00587   const MED_EN::medGeometryElement * types = mesh->getTypes( entity );
00588   for ( int t = 0; t < nbTypes; ++t )
00589   {
00590     // get connectivity
00591     if ( !mesh->existConnectivity( connType, entity ))
00592       continue;
00593     int nbCell = mesh->getNumberOfElements( entity, types[t] );
00594     const int *conn, *index;
00595     conn  = mesh->getConnectivity(connType, entity, types[t]);
00596     index = mesh->getConnectivityIndex(connType, entity);
00597     // find a cell containing the first of given nodes,
00598     // then check if the found cell contains all the given nodes
00599     const int firstNode = *nodes.begin();
00600     for ( int i = 0; i < nbCell; ++i )
00601     {
00602       for ( int j = index[i]-1; j < index[i+1]-1; ++j )
00603         if ( conn[j] == firstNode )
00604         {
00605           unsigned nbSame = 0;
00606           for ( j = index[i]-1; j < index[i+1]-1; ++j )
00607             nbSame += nodes.count( conn[j] );
00608           if ( nbSame == nodes.size() )
00609             return true;
00610           break;
00611         }
00612     }
00613   }
00614   return false;
00615 }
00616 
00618 void ParallelTopology::createFaceMapping(const MESHCollection& initial_collection,
00619                                          const MESHCollection& new_collection)
00620   //                     map<MED_EN::medGeometryElement,int*>& type_list,
00621   //                     map<MED_EN::medGeometryElement,int>& present_type_numbers,
00622   //                     int idomain
00623 
00624 {
00625   // containers for the new topology
00626   vector<int> new_counts(m_nb_domain,0);
00627   vector<int> domain_counts(m_nb_domain,0);
00628   const Topology* old_topology=initial_collection.getTopology();
00629   int nb_domain_old=old_topology->nbDomain();
00630   int global_index=old_topology->getFaceNumber();
00631   //cout << "nb faces = " << global_index << endl;
00632   set <pair<int, pair<int,int> > > global_treated;
00633 
00634   //definition of the d-1 constituent for the considered mesh dimension
00635   MED_EN::medEntityMesh constituent_entity;
00636   switch (m_mesh_dimension)
00637   {
00638   case 3:
00639     constituent_entity= MED_EN::MED_FACE;
00640     break;
00641   case 2:
00642     constituent_entity = MED_EN::MED_EDGE;
00643     break;
00644   }
00645 
00646   for (int iold=0; iold<nb_domain_old;iold++)
00647   {
00648     if ( !initial_collection.getMesh(iold) ) continue;
00649     int nbtotalface = initial_collection.getMesh(iold)->getNumberOfElements(constituent_entity,MED_EN::MED_ALL_ELEMENTS);
00650     SCRUTE_MED(nbtotalface);
00651     const int* face_conn = 0;
00652     const int* face_offset = 0;
00653     if (nbtotalface >0)
00654     {
00655       face_conn = initial_collection.getMesh(iold)->getConnectivity(
00656                                                                     MED_EN::MED_NODAL,constituent_entity,MED_EN::MED_ALL_ELEMENTS);
00657       face_offset = initial_collection.getMesh(iold)->getConnectivityIndex(MED_EN::MED_NODAL,constituent_entity);
00658     }
00659     for (int iface=0;iface<nbtotalface; iface++)
00660     {
00661       int global_face_number = old_topology->convertFaceToGlobal(iold,iface+1);
00662 
00663       //      int inode = face_offset[iface];
00664       for (int i=0; i<m_nb_domain; i++) domain_counts[i]=0;
00665       set <int> nodes;
00666       int nbnodes;
00667       {
00668         nbnodes=face_offset[iface+1]-face_offset[iface];
00669         for (int inode= face_offset[iface];inode < face_offset[iface+1]; inode++)
00670         {
00671           int node=face_conn[inode-1];
00672 
00673           int global = old_topology->convertNodeToGlobal(iold,node);
00674           //        cout << "global node "<<global<<"ip "<<iold<< "noeud"<<node<<endl;
00675           nodes.insert(global);
00676           typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00677           pair<mmiter,mmiter> range=m_node_glob_to_loc.equal_range(global);
00678 
00679           int ip;
00680           for (mmiter it=range.first; it !=range.second; it++)
00681           { 
00682             ip=(it->second).first;
00683             domain_counts[ip]++;
00684           }
00685         }
00686       }
00687 
00688       set<int>::const_iterator iter_node = nodes.begin();
00689       int numbers[3];
00690       numbers[2]=0; // for segments
00691       for (int i=0; i<nbnodes; i++)
00692       {
00693         numbers[i]=*iter_node;
00694         iter_node++;
00695       }
00696       set <pair<int, pair<int,int> > > ::iterator iter_triplets;
00697       pair<int, pair<int,int> > triplet = make_pair(numbers[0],make_pair(numbers[1],numbers[2]));
00698       iter_triplets=global_treated.find(triplet);
00699       if (iter_triplets==global_treated.end())
00700       {
00701         global_treated.insert(triplet);
00702         //  int nbnodes=face_offset[iface+1]-face_offset[iface];
00703         if (global_face_number == -1) 
00704         {
00705           global_index++;
00706           global_face_number=global_index;
00707 
00708         }
00709         //  SCRUTE_MED(nbnodes);
00710 
00711         for (int inew=0;inew<m_nb_domain;inew++)
00712         {
00713           //     SCRUTE_MED(domain_counts[inew]);
00714           if(domain_counts[inew]==nbnodes)
00715           {
00716             if ( !hasCellWithNodes( new_collection, inew, nodes ))
00717               continue; // 0020861: EDF 1387 MED: Result of medsplitter gives standalone triangles
00718 
00719             new_counts[inew]++;
00720             m_face_glob_to_loc.insert(make_pair(global_face_number,make_pair(inew,new_counts[inew])));
00721             //m_face_loc_to_glob.insert(make_pair(make_pair(inew,new_counts[inew]),global_face_number));
00722             m_face_loc_to_glob[inew].push_back(global_face_number);
00723           }
00724         }
00725         //cout << "global_face_number = " << global_face_number << endl;
00726       }
00727     }
00728   }
00729 
00730   for (int inew=0;inew<m_nb_domain;inew++)
00731   {
00732     m_nb_faces[inew]=new_counts[inew];
00733     MESSAGE_MED(" Nb faces ["<<inew<<"]="<<m_nb_faces[inew]);
00734   }
00735   MESSAGE_MED(" total number of faces"<<getFaceNumber());
00736 }
00737 
00739 void ParallelTopology::createFaceMapping2ndversion(const MESHCollection& initial_collection)
00740 
00741 {
00742   // containers for the new topology
00743   vector<int> new_counts(m_nb_domain,0);
00744   vector<int> domain_counts(m_nb_domain,0);
00745   const Topology* old_topology=initial_collection.getTopology();
00746   int nb_domain_old=old_topology->nbDomain();
00747   //int global_index=old_topology->getFaceNumber();
00748   //  set <pair<int, pair<int,int> > > global_treated;
00749 
00750   //definition of the d-1 constituent for the considered mesh dimension
00751   MED_EN::medEntityMesh constituent_entity;
00752   switch (m_mesh_dimension)
00753   {
00754   case 3:
00755     constituent_entity= MED_EN::MED_FACE;
00756     break;
00757   case 2:
00758     constituent_entity = MED_EN::MED_EDGE;
00759     break;
00760   }
00761 
00762   for (int iold=0; iold<nb_domain_old;iold++)
00763   {
00764     int nbcell=old_topology->getCellNumber(iold);
00765 
00766     const int* face_conn = initial_collection.getMesh(iold)->
00767       getConnectivity(MED_EN::MED_DESCENDING,MED_EN::MED_CELL,MED_EN::MED_ALL_ELEMENTS);
00768     const int* face_offset = initial_collection.getMesh(iold)->
00769       getConnectivityIndex(MED_EN::MED_DESCENDING,MED_EN::MED_CELL);
00770     MESSAGE_MED("end of connectivity calculation");
00771     set<int> global_treated;
00772     for (int icell=0; icell<nbcell; icell++)
00773     {
00774       int global_cell_number=old_topology->convertCellToGlobal(iold,icell+1);
00775       int inew = getCellDomainNumber(global_cell_number);
00776 
00777       for (int iface = face_offset[icell]; iface < face_offset[icell+1]; iface++)
00778       {
00779         int global_face_number=old_topology->convertFaceToGlobal(iold,abs(face_conn[iface-1]));
00780         if (global_treated.find (global_face_number)==global_treated.end())
00781         {
00782           new_counts[inew]++;
00783           m_face_glob_to_loc.insert(make_pair(global_face_number,make_pair(inew,new_counts[inew])));
00784 
00785           //m_face_loc_to_glob.insert(make_pair(make_pair(inew,new_counts[inew]),global_face_number));
00786           m_face_loc_to_glob[inew].push_back(global_face_number);
00787           global_treated.insert(global_face_number);
00788 
00789         }
00790       }
00791     }
00792   }
00793 
00794 
00795   for (int inew=0;inew<m_nb_domain;inew++)
00796   {
00797     m_nb_faces[inew]=new_counts[inew];
00798     //  cout << " Nb faces ["<<inew<<"]="<<m_nb_faces[inew]<<endl;
00799   }
00800   MESSAGE_MED(" total number of faces"<<getFaceNumber());
00801 }
00802 
00803 
00804 //replacing a table of global numbering with a table with local numberings
00805 // type_connectivity contains global connectivity for each type in input
00806 // type_connectivity contains local connectivity for each type in output
00807 void ParallelTopology::convertToLocal(map<MED_EN::medGeometryElement,int*>& type_connectivity,
00808                                       map<MED_EN::medGeometryElement,int>& present_type_numbers,
00809                                       int idomain,
00810                                       MED_EN::medEntityMesh entity)
00811 {
00812   int dimension;
00813   switch (entity)
00814   {
00815   case MED_EN::MED_CELL:
00816     dimension=m_mesh_dimension;
00817     break;
00818   case MED_EN::MED_FACE:
00819     dimension=2;
00820     break;
00821   case MED_EN::MED_EDGE:
00822     dimension=1;
00823     break;
00824   } 
00825 
00826   MED_EN::MESH_ENTITIES::const_iterator currentEntity;
00827   list<MED_EN::medGeometryElement>::const_iterator iter;
00828   currentEntity  = MED_EN::meshEntities.find(MED_EN::MED_CELL);
00829 
00830   for (iter = (*currentEntity).second.begin();iter != (*currentEntity).second.end(); iter++)
00831   {
00832     MED_EN::medGeometryElement type = (*iter);
00833     if (type/100 != dimension) continue;
00834     for (int inode=0; inode<present_type_numbers[type]*(type%100); inode++)
00835     {
00836       //      cout <<" inode :"<<inode<< " global = "<<type_connectivity[type][inode];
00837       int global = type_connectivity[type][inode];
00838       typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00839       pair<mmiter,mmiter> range=m_node_glob_to_loc.equal_range(global);
00840       for (mmiter it=range.first; it !=range.second; it++)
00841       {
00842         if ((it->second).first==idomain)
00843           type_connectivity[type][inode]=(it->second).second;
00844       }
00845       //      cout << " local = " <<type_connectivity[type][inode]<<endl;
00846     }
00847 
00848   }
00849 }
00850 
00851 //replacing a table of global numbering with a table with local numberings
00852 // type_connectivity contains global connectivity for each type in input
00853 // type_connectivity contains local connectivity for each type in output
00854 void ParallelTopology::convertToLocal2ndVersion(int* nodes, int nbnodes, int idomain)
00855 {
00856   for (int inode=0; inode<nbnodes; inode++)
00857   {
00858     //      cout <<" inode :"<<inode<< " global = "<<type_connectivity[type][inode];
00859     int global = nodes[inode];
00860     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::iterator mmiter;
00861     pair<mmiter,mmiter> range=m_node_glob_to_loc.equal_range(global);
00862     for (mmiter it=range.first; it !=range.second; it++)
00863     {
00864       if ((it->second).first==idomain)
00865         nodes[inode]=(it->second).second;
00866     }        
00867   }
00868 }
00869 
00870 
00872 void ParallelTopology::computeNodeNodeCorrespondencies(int idomain, vector<MEDMEM::MEDSKYLINEARRAY*>& corr) const
00873 {
00874   vector <int> sizes (m_nb_domain,0);
00875   vector <int*> values (m_nb_domain);
00876   for (int i=0; i<m_nb_domain; i++)
00877   {
00878     values[i]=0;
00879   }
00880 
00881   for (int inode=0; inode<m_nb_nodes[idomain]; inode++)
00882   {
00883     //int global = (m_node_loc_to_glob.find(make_pair(idomain,inode+1)))->second;
00884     int global = m_node_loc_to_glob[idomain][inode];
00885     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::const_iterator mm;
00886     pair<mm,mm> range = m_node_glob_to_loc.equal_range(global);
00887     for (mm it=range.first; it !=range.second; it++)
00888     {
00889       int id = (it->second).first;
00890       if (id !=idomain)
00891       {
00892         sizes[id]++;
00893       }
00894     }
00895   }
00896 
00897   for (int ip=0; ip<m_nb_domain; ip++)
00898   {
00899     if (sizes[ip]>0)
00900       values[ip]=new int[2*sizes[ip]];
00901     sizes[ip]=0;
00902   }
00903 
00904 
00905   for (int inode=0; inode<m_nb_nodes[idomain]; inode++)
00906   {
00907     //int global = (m_node_loc_to_glob.find(make_pair(idomain,inode+1)))->second;
00908     int global = m_node_loc_to_glob[idomain][inode];
00909     typedef INTERP_KERNEL::HashMultiMap<int,pair<int,int> >::const_iterator mm;
00910     pair<mm,mm> range = m_node_glob_to_loc.equal_range(global);
00911     for (mm it=range.first; it !=range.second; it++)
00912     {
00913       int id = (it->second).first;
00914       if (id !=idomain)
00915       {
00916         values[id][sizes[id]*2]=inode+1;
00917         values[id][sizes[id]*2+1]=(it->second).second;
00918         sizes[id]++;
00919       }
00920     }
00921   }
00922 
00923   for (int i=0; i< m_nb_domain; i++)
00924   {
00925     if (sizes[i]!=0)
00926     {
00927       int* index = new int[sizes[i]+1];
00928       for (int j=0; j<sizes[i]+1; j++)
00929         index[j]=j+1;
00930       corr[i]=new MEDMEM::MEDSKYLINEARRAY(sizes[i],2*sizes[i],index,values[i]);
00931       delete[] index;
00932       delete[] values[i];
00933     }
00934   }
00935 }
00936 
00937 //================================================================================
00946 //================================================================================
00947 
00948 void ParallelTopology::computeCellCellCorrespondencies(int idomain, vector<MEDMEM::MEDSKYLINEARRAY*>& corr, const Graph* graph) const
00949 {
00950   vector <int> sizes (m_nb_domain,0);
00951   vector <int*> values (m_nb_domain);
00952   for (int i=0; i<m_nb_domain; i++)
00953   {
00954     values[i]=0;
00955   }
00956 
00957   vector <INTERP_KERNEL::HashMultiMap<int,int> > cell_corresp;
00958   //TODO : remplacer int* par une map <int,int>
00959   //  vector<int*  > number_of_connections(m_nb_domain);
00960   //  vector<map<int,int> > number_of_connections;
00961   vector<map<int,int> > number_of_connections;
00962   cell_corresp.resize(m_nb_domain);
00963   number_of_connections.resize(m_nb_domain);
00964   //  for (int i=0; i<m_nb_domain; i++)
00965   //    {
00966   //      //    cell_corresp[i]=new multimap<int,int>;
00967   //      if (m_nb_cells[i] >0)
00968   //        {
00969   //          number_of_connections[i]=new int[m_nb_cells[idomain]];
00970   //          for (int j=0; j<m_nb_cells[idomain]; j++)
00971   //            number_of_connections[i][j]=0;
00972   //        }
00973   //    }
00974 
00975   const MEDMEM::MEDSKYLINEARRAY* skylinegraph = graph->getGraph();
00976 
00977   const int* index=skylinegraph->getIndex();
00978   const int* value=skylinegraph->getValue();
00979 
00980   TGlob2DomainLoc::const_iterator gl_do_lo_end = m_glob_to_loc.end();
00981 
00982   for (int icell=0; icell<m_nb_cells[idomain]; icell++)
00983   {
00984     int global= m_loc_to_glob[idomain][icell];
00985     for (int ii=index[global-1]-1; ii<index[global]-1; ii++)
00986     {
00987       int distant_global=value[ii];
00988 
00989       const pair<int,int>& local = m_glob_to_loc.find(distant_global)->second;
00990 
00991       if (local.first != idomain)
00992       {
00993         cell_corresp[local.first].insert(make_pair(icell+1,local.second));
00994         //              number_of_connections[local.first][icell]++;
00995         if (number_of_connections[local.first].find(icell)==number_of_connections[local.first].end())
00996           number_of_connections[local.first].insert(make_pair(icell,1));
00997         else
00998           number_of_connections[local.first][icell]++;
00999 
01000       }
01001     }
01002   }
01003 
01004   for (int inew=0; inew<m_nb_domain; inew++)
01005   {
01006     if (inew==idomain || number_of_connections[inew].empty()) continue;
01007 
01008     int* new_index=new int[m_nb_cells[idomain]+1];
01009     new_index[0]=1;
01010     for (int i=0; i<m_nb_cells[idomain]; i++)
01011     {
01012 
01013       if (number_of_connections[inew].find(i)!=number_of_connections[inew].end())
01014         new_index[i+1]=new_index[i]+number_of_connections[inew][i];
01015       else
01016         new_index[i+1]=new_index[i];
01017     }
01018     int* new_value;
01019     if (new_index[m_nb_cells[idomain]]-1 > 0)
01020       new_value=new int[new_index[m_nb_cells[idomain]]-1];
01021     else 
01022       new_value=0;
01023 
01024     int value_i=0;
01025 
01026     //                      INTERP_KERNEL::HashMultiMap<int,int>::iterator iter=cell_corresp[inew].begin();
01027 
01028     for (int i=0; i<m_nb_cells[idomain]; i++)
01029     {
01030       //          for (int j=new_index[i];j<new_index[i+1];j++,value_i++,iter++)
01031       //            new_value[value_i]=iter->second;
01032 
01033       typedef INTERP_KERNEL::HashMultiMap<int,int>::iterator mmiter;
01034       pair<mmiter,mmiter> range=cell_corresp[inew].equal_range(i+1);
01035       for (mmiter it=range.first; it!=range.second; it++)
01036       {
01037         new_value[value_i]=it->second;
01038         value_i++;
01039       }
01040     }
01041     if (value_i>0)    
01042     {
01043       corr[inew] = new MEDMEM::MEDSKYLINEARRAY(m_nb_cells[idomain], new_index[m_nb_cells[idomain]]-1, new_index,new_value,true);
01044     }
01045     else 
01046     {
01047       corr[inew]=0;
01048       if (new_value !=0) delete[]new_value;
01049       delete[]new_index;
01050     }   
01051 
01052   }
01053 
01054   //    for (int inew=0; inew<m_nb_domain; inew++)
01055   //    if (m_nb_cells[inew]>0)
01056   //      delete[] number_of_connections[inew];
01057 
01058 }
01059 
01060 //================================================================================
01064 //================================================================================
01065 
01066 int ParallelTopology::getMaxGlobalFace() const
01067 {
01068   int max = 0;
01069   TGlob2LocsMap::const_iterator g_l_l = m_face_glob_to_loc.begin();
01070   for ( ; g_l_l != m_face_glob_to_loc.end(); ++g_l_l )
01071     if ( g_l_l->first > max )
01072       max = g_l_l->first;
01073   return max;
01074 }
01075 
01076 void ParallelTopology::recreateFaceMapping(const TGeom2FacesByDomian& face_map)
01077 {
01078   m_face_glob_to_loc.clear();
01079   for (int i=0; i<m_nb_domain;i++)
01080     m_face_loc_to_glob[i].clear();
01081 
01082   for (int idomain=0; idomain<m_nb_domain; idomain++)
01083   {
01084     int ilocal=1;
01085     TGeom2Faces::const_iterator iter = face_map[idomain].begin();
01086     for (; iter != face_map[idomain].end(); iter++)
01087     {
01088       for (unsigned i=0; i< (iter->second).size(); i++)
01089       {
01090         MEDSPLITTER_FaceModel* face = (iter->second)[i];
01091         //cout << "global :"<<face->getGlobal()<<" - "<<ilocal<<endl;
01092         m_face_glob_to_loc.insert(make_pair(face->getGlobal(),make_pair(idomain,ilocal)));
01093         m_face_loc_to_glob[idomain].push_back(face->getGlobal());
01094         ilocal++;
01095       }
01096     }
01097     m_nb_faces[idomain] =ilocal-1;
01098   }
01099 }
01100 
01101 //================================================================================
01105 //================================================================================
01106 
01107 void ParallelTopology::recreateMappingAfterFusion(const vector<MEDMEM::MESH*>& meshes)
01108 {
01109   const char* LOC = "ParallelTopology::recreateMappingAfterFusion(): ";
01110 
01111   m_glob_to_loc.clear();
01112   m_node_glob_to_loc.clear();
01113   m_face_glob_to_loc.clear();
01114 
01115   for (int idomain=0; idomain<m_nb_domain; idomain++)
01116   {
01117     m_nb_cells[idomain] = m_nb_nodes[idomain] = m_nb_faces[idomain] = 0;
01118     m_loc_to_glob[idomain].clear();
01119     m_node_loc_to_glob[idomain].clear();
01120     m_face_loc_to_glob[idomain].clear();
01121     
01122     if ( !meshes[idomain]->getCoordinateptr() ) continue; // empty domian
01123 
01124     //creating cell maps
01125 
01126     m_nb_cells[idomain]=meshes[idomain]->getNumberOfElements(MED_EN::MED_CELL,
01127                                                              MED_EN::MED_ALL_ELEMENTS);
01128     if ( m_cell_loc_to_glob_fuse[idomain].size() != m_nb_cells[idomain] )
01129       throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused cells");
01130 
01131     m_loc_to_glob[idomain].swap( m_cell_loc_to_glob_fuse[idomain] );
01132 
01133     for (int i=0; i<m_nb_cells[idomain]; i++)
01134     {
01135       int global=m_loc_to_glob[idomain][i];
01136       m_glob_to_loc[global]=make_pair(idomain,i+1);
01137     }
01138 
01139     //creating node maps
01140 
01141     m_nb_nodes[idomain]=meshes[idomain]->getNumberOfNodes();
01142     m_node_loc_to_glob[idomain] = ((MEDMEM::MeshFuse*)meshes[idomain])->getNodeNumbers();
01143     if ( m_node_loc_to_glob[idomain].size() != m_nb_nodes[idomain] )
01144       throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused nodes");
01145 
01146     // setting mappings for all nodes
01147     for (int inode=0; inode<m_nb_nodes[idomain]; inode++)
01148     {
01149       int global_number=m_node_loc_to_glob[idomain][inode];
01150       m_node_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,inode+1)));
01151     }
01152 
01153 
01154     //creating dimension d-1 component mappings
01155 
01156     MED_EN::medEntityMesh constituent_entity =
01157       (m_mesh_dimension == 3 ? MED_EN::MED_FACE : MED_EN::MED_EDGE );
01158     m_nb_faces[idomain] = meshes[idomain]->getNumberOfElements(constituent_entity,
01159                                                                MED_EN::MED_ALL_ELEMENTS);
01160     if ( m_face_loc_to_glob_fuse[idomain].size() != m_nb_faces[idomain] )
01161       throw MED_EXCEPTION(MEDMEM::STRING(LOC)<<" invalid nb of fused faces of domain "<< idomain
01162                           << ": expect " << m_nb_faces[idomain]
01163                           << " but have " << m_face_loc_to_glob_fuse[idomain].size());
01164 
01165     m_face_loc_to_glob[idomain].swap( m_face_loc_to_glob_fuse[idomain] );
01166 
01167     for (int iface=0; iface<m_nb_faces[idomain]; iface++)
01168     {
01169       int global_number=m_face_loc_to_glob[idomain][iface];
01170       m_face_glob_to_loc.insert(make_pair(global_number,make_pair(idomain,iface+1)));
01171     }
01172   }
01173 
01174 }