Salome HOME
Merge from tetra_hpc branch
[plugins/ghs3dprlplugin.git] / src / tepal2med / ghs3dprl_mesh_wrap.cxx
index 1ed2e55aa3211530f5eab7b4e2c42b3c50014f91..8f557f4fad710ed7d56c082996e563e6fb648e4c 100755 (executable)
@@ -448,6 +448,40 @@ bool CVW_FindString(const std::string &str,std::fstream &Ff, long &count)
    return true;
 }
 
+//************************************
+bool CVW_FindStringInFirstLines(const std::string &str,std::fstream &Ff, long &maxline)
+//find in file maximum maxline first lines with string str in first position of line
+//converts count value expected after " " in line found
+{
+   std::string line;
+   long nbLine=0;
+   do
+   {
+      if (getline(Ff,line))
+      {
+         nbLine++;
+         if (line[0]==str[0]) //faster
+         {
+            if (line.find(str)==0)
+            {
+            return true;
+            }
+         }
+         if (nbLine>=maxline)
+         {
+            std::cerr<<"Problem line '"<<str<<"' not found in first "<<maxline<<" lines of file\n"<<std::endl;
+            return false;
+         }
+      }
+      else
+      {
+         std::cerr<<"Problem line '"<<str<<"' not found in all file\n"<<std::endl;
+         return false;
+      }
+   } while (1);
+   return true;
+}
+
 //************************************
 bool ghs3dprl_mesh_wrap::ReadFileMSGnew(const QString FileName)
 //read file .glo with no parser xml because big file (volume)
@@ -519,6 +553,34 @@ bool ghs3dprl_mesh_wrap::ReadFileMSGnew(const QString FileName)
    return true;
 }
 
+//************************************
+bool ghs3dprl_mesh_wrap::TestExistingFileMESHnew(const QString FileName)
+//read file .glo with no parser xml because big file (volume)
+//no read of <receive> for speed (and so no test)
+{
+   QString tmp;
+   std::fstream Ff(FileName.toLatin1().constData(),std::ios_base::in);
+   std::string line;
+   long i,count,meshversion,maxline;
+   bool ok;
+
+   if (!Ff.is_open())
+   {
+      std::cerr<<"Problem File '"<<FileName.toLatin1().constData()<<"' not open\n"<<std::endl;
+      return false;
+   }
+
+   //Lit les donn�s au debut du fichier, 1 lignes maxi:
+   maxline=1;
+   if (!CVW_FindStringInFirstLines("MeshVersionFormatted",Ff,maxline)) return false;
+   if (verbose>2) std::cout<<"MeshVersionFormatted_"<<this->nofile<<" ok"<<std::endl;
+
+   //Ferme le fichier :
+   Ff.close();
+   this->nbfiles++;
+   return true;
+}
+
 ///************************************
 bool ghs3dprl_mesh_wrap::ReadFileGLO(const QString FileName)
 //read file .glo with no parser xml because big file (volume)
@@ -645,6 +707,171 @@ bool ghs3dprl_mesh_wrap::ReadFileFACES(const QString FileName)
    return true;
 }
 
+//************************************
+bool ghs3dprl_mesh_wrap::ReadFileMESH(const QString FileName)
+//read file .mesh from tetra_hpc (with volume)
+{
+   std::string line("                                            ");
+   QString tmp;
+   std::fstream Ff(FileName.toLatin1().constData(),std::ios_base::in);
+   long Mversion=0,Mdim=0,Mvert=0,Mtria=0,Mtetra=0;
+   med_int garbage;
+   //long ne,np,npfixe,subnumber,reste;
+   bool ok;
+   if (verbose>=6)std::cout<<"Read File '"<<FileName.toLatin1().constData()<<std::endl;
+   if (!Ff.is_open()){
+      std::cerr<<"Problem File '"<<FileName.toLatin1().constData()<<"' not open\n"<<std::endl;
+      return false;
+   }
+
+   //lit les données :
+   if (getline(Ff,line) && (line.find("MeshVersionFormatted")==0))
+   {
+      tmp=line.c_str();
+      Mversion=tmp.section(' ',1,1).toLong(&ok);
+   }
+   else
+   {
+      std::cerr<<"Problem on line 1 of file: 'MeshVersionFormatted' expected"<<std::endl;
+      return false;
+   }
+   
+   getline(Ff,line);
+
+   if (getline(Ff,line) && (line.find("Dimension 3")==0))
+   {
+      tmp=line.c_str();
+      Mdim=tmp.section(' ',1,1).toLong(&ok);
+   }
+   else
+   {
+      std::cerr<<"Problem on line 3 of file: Dimension 3 expected"<<std::endl;
+      return false;
+   }
+   
+   getline(Ff,line);
+
+   if (!(getline(Ff,line) && (line.find("Vertices")==0)))
+   {
+      std::cerr<<"Problem on line 5 of file: 'Vertices' expected"<<std::endl;
+      return false;
+   }
+   if (getline(Ff,line))
+   {
+      tmp=line.c_str();
+      Mvert=tmp.toLong(&ok);
+   }
+   else
+   {
+      std::cerr<<"Problem on line 6 of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+   if (Mvert<=0)
+   {
+      std::cerr<<"Problem on line 6 of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+
+   med_float *tmflo=new med_float[Mvert*3];
+   for (int i=0; i<Mvert*3; i=i+3) Ff>>tmflo[i]>>tmflo[i+1]>>tmflo[i+2]>>garbage;
+   if (verbose>4) std::cout<<"Vertices "<<tmflo[0]<<" "<<tmflo[1]<<"... "<<tmflo[Mvert*3-1]<<std::endl;
+
+   CVWtab *montab=new CVWtab(Mvert*3,tmflo);
+   tmp=tmp.sprintf("NB%ld VC",this->nofile);
+   ok=this->insert_key(tmp,montab);
+
+   getline(Ff,line);
+   getline(Ff,line);
+
+   if (!(getline(Ff,line) && (line.find("Triangles")==0)))
+   {
+      std::cerr<<"Problem on line 'Triangles' of file: not found"<<line<<std::endl;
+      return false;
+   }
+   if (getline(Ff,line))
+   {
+      tmp=line.c_str();
+      Mtria=tmp.toLong(&ok);
+   }
+   else
+   {
+      std::cerr<<"Problem on line 'Triangles' of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+   if (Mtria<=0)
+   {
+      std::cerr<<"Problem on line 'Triangles' of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+
+   med_int *tmint=new med_int[Mtria*7]; //*7 as older files .faces
+   for (int i=0; i<Mtria*7; i=i+7) {
+     Ff>>tmint[i]>>tmint[i+1]>>tmint[i+2]>>garbage;
+     tmint[i+3]=0;
+     tmint[i+4]=0;
+     tmint[i+5]=0;
+     tmint[i+6]=0;
+   }
+   if (verbose>4) std::cout<<"Triangles "<<tmint[0]<<" "<<tmint[1]<<"... "<<
+       tmint[Mtria*7-5]<<" "<<tmint[Mtria*7-4]<<" "<<tmint[Mtria*7-3]<<" "<<
+      tmint[Mtria*7-2]<<" "<<tmint[Mtria*7-1]<<std::endl;
+
+   montab=new CVWtab(Mtria*7,tmint);
+   tmp=tmp.sprintf("FC%ld",this->nofile);
+   ok=this->insert_key(tmp,montab);
+   
+   getline(Ff,line);
+   getline(Ff,line);
+
+   if (!(getline(Ff,line) && (line.find("Tetrahedra")==0)))
+   {
+      std::cerr<<"Problem on line 'Tetrahedra' of file: not found"<<std::endl;
+      return false;
+   }
+   if (getline(Ff,line))
+   {
+      tmp=line.c_str();
+      Mtetra=tmp.toLong(&ok);
+   }
+   else
+   {
+      std::cerr<<"Problem on line 'Tetrahedra' of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+   if (Mtetra<=0)
+   {
+      std::cerr<<"Problem on line 'Tetrahedra' of file: a positive integer is expected"<<std::endl;
+      return false;
+   }
+
+   if (verbose>=2)
+   {
+      std::cout<<"MeshVersionFormatted="<<Mversion<<std::endl;
+      std::cout<<"MeshDimension="<<Mdim<<std::endl;
+      std::cout<<"MeshNumberOfVertices="<<Mvert<<std::endl;
+      std::cout<<"MeshNumberOfTriangles="<<Mtria<<std::endl;
+      std::cout<<"MeshNumberOfTetrahedra="<<Mtetra<<std::endl;
+   }
+
+   tmint=new med_int[Mtetra*4];
+   for (int i=0; i<Mtetra*4; i=i+4) Ff>>tmint[i]>>tmint[i+1]>>tmint[i+2]>>tmint[i+3]>>garbage;
+   if (verbose>4) std::cout<<"Tetrahedra "<<tmint[0]<<" "<<tmint[1]<<"... "<<tmint[Mtetra*4-1]<<std::endl;
+
+   montab=new CVWtab(Mtetra*4,tmint);
+   tmp=tmp.sprintf("NB%ld EV",this->nofile);
+   ok=this->insert_key(tmp,montab);
+
+   //swap on file if too big for memory in one cpu
+   //default 1GOctet/8(for double)/10(for arrays in memory at the same time)
+   if (Mvert*3>this->nbelem_limit_swap)
+     this->SwapOutOfMemory_key_mesh_wrap(QRegExp("NB",Qt::CaseSensitive,QRegExp::RegExp));
+   //beware record 6 lenght 1
+   //ferme le fichier :
+   Ff.close();
+   this->nbfiles++;
+   return true;
+}
+
 //************************************
 bool ghs3dprl_mesh_wrap::ReadFileNOBOITE(const QString FileName)
 //read file .noboite (volume)
@@ -1370,7 +1597,7 @@ bool ghs3dprl_mesh_wrap::Write_masterxmlMEDfile()
    //Creating child nodes
    //Version tag
    med_int majeur,mineur,release;
-   //Quelle version de MED est utilis
+   //Quelle version de MED est utilisee
    MEDlibraryNumVersion(&majeur,&mineur,&release);
    if (verbose>0) fprintf(stdout,"Files write with MED V%d.%d.%d\n",majeur,mineur,release);
    node = xmlNewChild(root_node, 0, BAD_CAST "version",0);
@@ -1381,7 +1608,7 @@ bool ghs3dprl_mesh_wrap::Write_masterxmlMEDfile()
 
    //Description tag
    node = xmlNewChild(root_node,0, BAD_CAST "description",0);
-   xmlNewProp(node, BAD_CAST "what", BAD_CAST "tetrahedral mesh by tepal");
+   xmlNewProp(node, BAD_CAST "what", BAD_CAST "tetrahedral mesh by MeshGems-Tetra-hpc (formerly tepal)");
 #ifdef WIN32
   SYSTEMTIME  present;
   GetLocalTime ( &present );
@@ -1522,7 +1749,8 @@ bool ghs3dprl_mesh_wrap::Write_MEDfiles_v2(bool deletekeys)
                                                 cmd.toLatin1().constData()<<"> does not exist\n"; }
 
    //define family 0 if not existing, no groups
-   families.add("0","");
+   //La famille FAMILLE_ZERO n'a pas été trouvée, elle est obligatoire
+   families.add("0","FAMILLE_ZERO");
    //define family Group_of_New_Nodes (which not exists before tetrahedra)
    famallnodes=0;
    if (QString("All_Nodes").contains(deletegroups)==0){
@@ -1647,6 +1875,7 @@ bool ghs3dprl_mesh_wrap::Write_MEDfiles_v2(bool deletekeys)
       MEDfileClose(fid); //but loop on others domains
 
    }
+   
    MEDfileClose(fidjoint); //no error
    if (verbose>0)std::cout<<"\nTotalNumberOftetrahedra="<<nbtetrastotal<<std::endl;
 
@@ -1667,8 +1896,14 @@ bool ghs3dprl_mesh_wrap::idom_nodes()
       key=key.sprintf("NB%d VC",idom); //files.NoBoite Vertex Coordinates
       tab=this->restore_key(key); //tab1=this->mestab[key1];
       if (!tab) {
-         tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".noboite";
-         ok=this->ReadFileNOBOITE(tmp);
+         if (!for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".noboite";
+            ok=this->ReadFileNOBOITE(tmp);
+         }
+         if (for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format_tetra.toLatin1().constData(),idom)+".mesh";
+            ok=this->ReadFileMESH(tmp);
+         }
          tab=this->restore_key(key); //tab1=this->mestab[key1];
          if (!tab) return false;
       }
@@ -1690,8 +1925,14 @@ bool ghs3dprl_mesh_wrap::idom_nodes()
       key1=key1.sprintf("GL%d VE",idom); //global numerotation 
       tab1=this->restore_key(key1); //tab1=this->mestab[key1];
       if (!tab1) {
-         tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".glo";
-         ok=this->ReadFileGLO(tmp);
+         if (!for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".glo";
+            ok=this->ReadFileGLO(tmp);
+         }
+         if (for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format_tetra.toLatin1().constData(),idom)+".glo";
+            ok=this->ReadFileGLO(tmp);
+         }
          if (!ok) {std::cerr<<"Problem file "<<tmp.toLatin1().constData()<<std::endl; return false;}
          tab1=this->restore_key(key1); //tab1=this->mestab[key1];
          if (!tab1) return false;
@@ -1942,17 +2183,14 @@ bool ghs3dprl_mesh_wrap::idom_joints()
       //qt3 if (this->nb_key_mesh_wrap(QRegExp(tmp,true,true))<=0) {
       if (this->nb_key_mesh_wrap(QRegExp(tmp,Qt::CaseSensitive,QRegExp::RegExp))<=0) {
          this->nofile=idom;
-         /*old version with xml parser too slow
-         ghs3dprl_msg_parser handler;
-         handler.mailw=this;
-         QXmlSimpleReader reader;
-         reader.setContentHandler(&handler);
-         tmp=pathini+casename+tmp.sprintf(format,nbfilestot,idom)+".msg";
-         QFile File(tmp);
-         QXmlInputSource source(&File);
-         reader.parse(source);
-         File.close();*/
-         tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".msg";
+         
+         if (!for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,idom)+".msg";
+         }
+         if (for_tetrahpc) {
+            tmp=pathini+casename+tmp.sprintf(format_tetra.toLatin1().constData(),idom)+".msg";
+         }
+
          ok=this->ReadFileMSGnew(tmp);
          if (!ok) {
             std::cerr<<"Problem in file "<<tmp.toLatin1().constData()<<std::endl;
@@ -1973,7 +2211,14 @@ bool ghs3dprl_mesh_wrap::idom_joints()
          tab2=this->restore_key(key1);
          //if not yet loaded (first time) try to load necessary file msg of neighbourg
          if (!tab2) {
-            tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,ineig)+".msg";
+            
+            if (!for_tetrahpc) {
+               tmp=pathini+casename+tmp.sprintf(format.toLatin1().constData(),nbfilestot,ineig)+".msg";
+            }
+            if (for_tetrahpc) {
+               tmp=pathini+casename+tmp.sprintf(format_tetra.toLatin1().constData(),ineig)+".msg";
+            }
+            
             this->nofile=ineig; //neighbourg file
             ok=this->ReadFileMSGnew(tmp);
             this->nofile=idom;  //restaure initial domain
@@ -2019,14 +2264,14 @@ bool ghs3dprl_mesh_wrap::idom_joints()
          //two indices for one correspondence
          arrayi=new med_int[nb*2];
          arraynodes=new med_float[nbnodesneig*3];  //for file DOMAIN_join.med
-         inodes=new med_int[nbnodes];            //for file DOMAIN_join.med
+         inodes=new med_int[nbnodes];              //for file DOMAIN_join.med
          med_int * fammore=new med_int[nbnodes];
-         for (i=0; i<nbnodes ; i++) {fammore[i]=0; inodes[i]=-1;}
-         ii=0; jj=0; k=0;
+         for (i=0; i<nbnodes ; i++) {fammore[i]=0; inodes[i]=-2;}  //precautions
+         ii=0; jj=0; k=1;
          for (i=0; i<nb ; i++){
             //no need because <send> equals <receive> tab1->tmint[i]==tab2->tmint[i]
-            j=tab1->tmint[i]-1; //contents of tab1 1 to nb
-            inodes[j]=k; k++;   //indices 0->n-1 of nodes of joint from nodes of domain
+            j=tab1->tmint[i]-1; //contents of tab1 1->nb, j 0->nb-1
+            inodes[j]=k; k++;   //contents of inodes 1->n ,nodes of joint from nodes of domain
             arraynodes[jj]=tab->tmflo[j*3]; jj++;
             arraynodes[jj]=tab->tmflo[j*3+1]; jj++;
             arraynodes[jj]=tab->tmflo[j*3+2]; jj++;
@@ -2061,7 +2306,7 @@ bool ghs3dprl_mesh_wrap::idom_joints()
          key1=key1.sprintf("MS%d NE%d FA SE",ineig,idom); //SE or RE identicals
          tab2=this->restore_key(key1);
          if (!tab2) std::cerr<<"Problem existing triangles of joint in domain "<<idom<<
-                          " with none in neighbourg "<<ineig<<" files .msg"<<std::endl;
+                               " with none in neighbourg "<<ineig<<" files .msg"<<std::endl;
          nb=tab1->size; nbtria3neig=tab2->size;
          if (nb!=nbtria3neig) {
             std::cerr<<"Problem in file "<<tmp.toLatin1().constData()<<
@@ -2080,24 +2325,43 @@ bool ghs3dprl_mesh_wrap::idom_joints()
          tab=this->restore_key(key); //tab1=this->mestab[key1];
 
          med_int nb=tab1->size; nbtria3neig=nb;
-         if (verbose>4)
-            std::cout<<"NumberOfTrianglesOfJoint_"<<idom<<"_"<<ineig<<"="<<nb<<std::endl;
-         arrayi=new med_int[nb*2];
+         //if (verbose>=0) std::cout<<"NumberOfTrianglesOfJoint_"<<idom<<"_"<<ineig<<"="<<nb<<std::endl;
+         arrayi=new med_int[nb*2]; //correspondance indices triangles in 2 domains
          arrayfaces=new med_int[nbtria3neig*3];  //for file DOMAIN_join.med
+         for (i=0; i<nbtria3neig*3 ; i++) arrayfaces[i]=-1; //precaution
          fammore=new med_int[nbtria3];
          for (i=0; i<nbtria3 ; i++) fammore[i]=0;
          ii=0; jj=0;
          for (i=0; i<nb ; i++){
             arrayi[ii]=tab1->tmint[i]; ii++;
+            arrayi[ii]=tab2->tmint[i]; ii++; //correspondance indices triangles in 2 domains
+       
             fammore[tab1->tmint[i]-1]=famjoint;
             //famtria3[tab1->tmint[i]-1]=famjoint;
-            arrayi[ii]=tab2->tmint[i]; ii++;
-            //std::cout<<arrayi[ii-1]<<"="<<arrayi[ii-2]<<std::endl;
-            k=(tab1->tmint[i]-1)*7; //indice of node connectivity
-            arrayfaces[jj]=inodes[tab->tmint[k]-1]+1; jj++;
-            arrayfaces[jj]=inodes[tab->tmint[k+1]-1]+1; jj++;
-            arrayfaces[jj]=inodes[tab->tmint[k+2]-1]+1; jj++;
+            
+            k=tab1->tmint[i]-1; //indice of node connectivity
+            //std::cout<<"k="<<k<<std::endl;
+            k=k*7; //indice of node connectivity in tab of triangles
+            
+            arrayfaces[jj]=inodes[tab->tmint[k]-1]; jj++;
+            arrayfaces[jj]=inodes[tab->tmint[k+1]-1]; jj++;
+            arrayfaces[jj]=inodes[tab->tmint[k+2]-1]; jj++;
+         }
+         int happens=0;
+         for (i=0; i<nbtria3neig*3 ; i++) {
+           if (arrayfaces[i]<=0) {
+             std::cerr<<"Problem file X_joints.med unknown node in joint "<<idom<<"_"<<ineig<<" face "<<i/3+1<<std::endl; //precaution
+             happens=1;
+           }
+         }
+         /*TODO DEBUG may be bug distene?
+         if (happens==1) {
+            std::cout<<"\nNumberOfTrianglesOfJoint_"<<idom<<"_"<<ineig<<"="<<nb<<std::endl;
+            for (i=0; i<nbnodes ; i++) std::cout<<"inode i "<<i+1<<" "<<inodes[i]<<std::endl;
+            for (i=0; i<tab1->size ; i++) std::cout<<"triangle i "<<i+1<<" "<<tab1->tmint[i]<<std::endl;
+            for (i=0; i<tab->size ; i=i+7) std::cout<<"conn i "<<i/7+1<<" : "<<tab->tmint[i]<<" "<<tab->tmint[i+1]<<" "<<tab->tmint[i+2]<<std::endl;
          }
+         */
          if (namejoint.contains(deletegroups)==0){
             ok=set_one_more_family(famtria3,fammore,nbtria3);
          }
@@ -2126,7 +2390,10 @@ bool ghs3dprl_mesh_wrap::idom_joints()
           err=MEDmeshEntityNumberWr(fidjoint,namejnt,MED_NO_DT,MED_NO_IT,MED_NODE,MED_UNDEF_GEOMETRY_TYPE,nbnodesneig,arrayi);
           delete[] arrayi;
           if (err<0) std::cerr<<"Problem MEDmeshEntityNumberWr of nodes "<<namejnt<<std::endl;
-
+          //families zero in file fidjoint ???
+          //La famille FAMILLE_ZERO n'a pas été trouvée, elle est obligatoire
+          nb=create_family_zero(fidjoint,namejnt);
+          
           //write tria3
           if (nbtria3neig>0) {
            //for (i=0; i<nbtria3neig ; i++) std::cout<<i+1<<" "<<
@@ -2266,10 +2533,29 @@ med_int ghs3dprl_mesh_wrap::create_families(med_idt fid, int sign)
       }
       if (verbose>5)std::cout<<std::endl;
       err=MEDfamilyCr(fid,nomfinal,nomfam,num,ngro,gro);
-      if (err<0) std::cerr<<"Problem MEDfamilyCr"<<std::endl;
       delete[] gro;
       if (err<0) std::cerr<<"Problem MEDfamilyCr of "<<nomfam<<std::endl;
    }
    return ires;
 }
 
+med_int ghs3dprl_mesh_wrap::create_family_zero(med_idt fid, QString nameMesh)
+{
+   med_int pas,ires;
+   ires=0;
+   char nomfam[MED_NAME_SIZE+1]="FAMILLE_ZERO";  //it.current()->name;
+   char attdes[MED_COMMENT_SIZE+1]="_NO_DESCRIPTION";
+   
+   char *gro;
+   med_int i,attide=1,attval=1,natt=1,num=0,ngro=0;
+   
+   gro=new char[MED_LNAME_SIZE*ngro+2];
+   if (verbose>3)std::cout<<"\ncreate_family_ZERO "<<nameMesh.toLatin1().constData()<<std::endl;
+   err=MEDfamilyCr(fid,nameMesh.toLatin1().constData(),nomfam,num,ngro,gro);
+   if (err<0) std::cerr<<"Problem MEDfamilyCr FAMILLE_ZERO of "<<nameMesh.toLatin1().constData()<<std::endl;
+   delete[] gro;
+   return ires;
+}
+
+
+