When loading a netgen-mesh, take geometry-info from mesh-file over global geometry.

libsrc/interface/nginterface.cpp

@@ -157,7 +157,7 @@ void Ng_LoadMesh (const char * filename, ngcore::NgMPI_Comm comm)
     SetGlobalMesh (mesh);
 
     // make string from rest of file (for geometry info!)
-    if(!ng_geometry) {
+    // (this might be empty, in which case we take the global ng_geometry)
     stringstream geom_part;
     geom_part << infile->rdbuf();
     string geom_part_string = geom_part.str();
@@ -165,7 +165,7 @@ void Ng_LoadMesh (const char * filename, ngcore::NgMPI_Comm comm)
     // buf = new char[strs];
     buf.SetSize(strs);
     memcpy(&buf[0], geom_part_string.c_str(), strs*sizeof(char));
-    }
+
     delete infile;
 
     if (ntasks > 1)
@@ -238,32 +238,23 @@ void Ng_LoadMesh (const char * filename, ngcore::NgMPI_Comm comm)
     mesh->SendRecvMesh();
   }
 
-  if(!ng_geometry && ntasks>1) {
+  if(ntasks>1) {
 #ifdef PARALLEL
     /** Scatter the geometry-string (no dummy-implementation in mpi_interface) **/
+    int strs = buf.Size();
+    MyMPI_Bcast(strs, comm);
+    if(strs>0)
       MyMPI_Bcast(buf, comm);
 #endif
   }
 
-  if(!ng_geometry) {
+  shared_ptr<NetgenGeometry> geo;
-    infile = new istringstream(string((const char*)&buf[0], (size_t)strs));
+  if(buf.Size()) { // if we had geom-info in the file, take it
-    // delete[] buf;
+    istringstream geom_infile(string((const char*)&buf[0], buf.Size()));
-    for (int i = 0; i < geometryregister.Size(); i++)
+    geo = geometryregister.LoadFromMeshFile(geom_infile);
-      {
-	NetgenGeometry * hgeom = geometryregister[i]->LoadFromMeshFile (*infile);
-	if (hgeom)
-	  {
-	    ng_geometry.reset (hgeom);
-	    mesh->SetGeometry(ng_geometry);
-	    break;
   }
-      }
+  if(geo!=nullptr) mesh->SetGeometry(geo);
-  }
+  else if(ng_geometry!=nullptr) mesh->SetGeometry(ng_geometry);
-  /** Dummy Geometry if we still could not find any geometry info! **/
-  // if (!ng_geometry)
-  //   ng_geometry = make_shared<NetgenGeometry>();
-  if(ng_geometry)
-    mesh->SetGeometry(ng_geometry);
 }
 
 void Ng_LoadMeshFromString (const char * mesh_as_string)

libsrc/meshing/python_mesh.cpp

@@ -578,66 +578,139 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
         return mesh;
       })
     .def("Load",  FunctionPointer 
-	 ([](Mesh & self, const string & filename)
+	 ([](shared_ptr<Mesh> self, const string & filename)
 	  {
+
+	    auto comm = self->GetCommunicator();
+	    int id = comm.Rank();
+	    int ntasks = comm.Size();
+	    auto & mesh = self;
+
+	    {
+	      ifstream infile(filename.c_str());
+	      if(!infile.good())
+		throw NgException(string("Error opening file ") + filename);
+	    }
+
+	    if ( filename.find(".vol") == string::npos )
+	      {
+		if(ntasks>1)
+		  throw NgException("Not sure what to do with this?? Does this work with MPI??");
+		mesh->SetCommunicator(comm);
+		ReadFile(*mesh,filename.c_str());
+		//mesh->SetGlobalH (mparam.maxh);
+		//mesh->CalcLocalH();
+		return;
+	      }
+
 	    istream * infile;
+	    Array<char> buf; // for distributing geometry!
+	    int strs;
 
-	    NgMPI_Comm comm = self.GetCommunicator();
-	    id = comm.Rank();
-	    ntasks = comm.Size();
-
-#ifdef PARALLEL
-	    char* buf = nullptr;
-	    int strs = 0;
 	    if( id == 0) {
-#endif
+
-            if (filename.find(".vol.gz") != string::npos)
+	      if (filename.substr (filename.length()-3, 3) == ".gz")
 		infile = new igzstream (filename.c_str());
 	      else
 		infile = new ifstream (filename.c_str());
-	    // ifstream input(filename);
+	      mesh -> Load(*infile);
-#ifdef PARALLEL
-	    //still inside id==0-bracket...
-	        self.Load(*infile);	      
-		self.Distribute();
 
-		/** Copy the rest of the file into a string (for geometry) **/
+	      // make string from rest of file (for geometry info!)
+	      // (this might be empty, in which case we take the global ng_geometry)
 	      stringstream geom_part;
 	      geom_part << infile->rdbuf();
 	      string geom_part_string = geom_part.str();
 	      strs = geom_part_string.size();
-		buf = new char[strs];
+	      // buf = new char[strs];
-		memcpy(buf, geom_part_string.c_str(), strs*sizeof(char));
+	      buf.SetSize(strs);
+	      memcpy(&buf[0], geom_part_string.c_str(), strs*sizeof(char));
+
+	      delete infile;
+
+	      if (ntasks > 1)
+		{
+
+		  char * weightsfilename = new char [filename.size()+1];
+		  strcpy (weightsfilename, filename.c_str());            
+		  weightsfilename[strlen (weightsfilename)-3] = 'w';
+		  weightsfilename[strlen (weightsfilename)-2] = 'e';
+		  weightsfilename[strlen (weightsfilename)-1] = 'i';
+
+		  ifstream weightsfile(weightsfilename);      
+		  delete [] weightsfilename;  
+	  
+		  if (!(weightsfile.good()))
+		    {
+		      // cout << "regular distribute" << endl;
+		      mesh -> Distribute();
 		    }
+		  else
+		    {
+		      char str[20];   
+		      bool endfile = false;
+		      int n, dummy;
+	      
+		      Array<int> segment_weights;
+		      Array<int> surface_weights;
+		      Array<int> volume_weights;
+	      
+		      while (weightsfile.good() && !endfile)
+			{
+			  weightsfile >> str;
+		  
+			  if (strcmp (str, "edgeweights") == 0)
+			    {
+			      weightsfile >> n;
+			      segment_weights.SetSize(n);
+			      for (int i = 0; i < n; i++)
+				weightsfile >> dummy >> segment_weights[i];
+			    }
+		  
+			  if (strcmp (str, "surfaceweights") == 0)
+			    {
+			      weightsfile >> n;
+			      surface_weights.SetSize(n);
+			      for (int i=0; i<n; i++)
+				weightsfile >> dummy >> surface_weights[i];
+			    }
+		  
+			  if (strcmp (str, "volumeweights") == 0)
+			    {
+			      weightsfile >> n;
+			      volume_weights.SetSize(n);
+			      for (int i=0; i<n; i++)
+				weightsfile >> dummy >> volume_weights[i];
+			    }
+		  
+			  if (strcmp (str, "endfile") == 0)
+			    endfile = true;  
+			}     
+	      
+		      mesh -> Distribute(volume_weights, surface_weights, segment_weights);
+		    }
+		} // ntasks>1 end
+	    } // id==0 end
 	    else {
-	      self.SendRecvMesh();
+	      mesh->SendRecvMesh();
 	    }
 
-	    /** Scatter the geometry-string **/
+	    if(ntasks>1) {
-	    MPI_Bcast(&strs, 1, MPI_INT, 0, comm); 
+#ifdef PARALLEL
-	    if(id!=0)
+	      /** Scatter the geometry-string (no dummy-implementation in mpi_interface) **/
-	      buf = new char[strs];
+	      int strs = buf.Size();
-	    MPI_Bcast(buf, strs, MPI_CHAR, 0, comm);
+	      MyMPI_Bcast(strs, comm);
-	    if(id==0)
+	      if(strs>0)
-	      delete infile;
+		MyMPI_Bcast(buf, comm);
-	    infile = new istringstream(string((const char*)buf, (size_t)strs));
-	    delete[] buf;
-	    
-#else
-	    self.Load(*infile);
 #endif
-	    for (int i = 0; i < geometryregister.Size(); i++)
-	      {
-		NetgenGeometry * hgeom = geometryregister[i]->LoadFromMeshFile (*infile);
-		if (hgeom)
-		  {
-		    ng_geometry.reset (hgeom);
-                    self.SetGeometry(ng_geometry);
-		    break;
 	    }
+
+	    shared_ptr<NetgenGeometry> geo;
+	    if(buf.Size()) { // if we had geom-info in the file, take it
+	      istringstream geom_infile(string((const char*)&buf[0], buf.Size()));
+	      geo = geometryregister.LoadFromMeshFile(geom_infile);
 	    }
-	    self.SetGeometry(ng_geometry);
+	    if(geo!=nullptr) mesh->SetGeometry(geo);
-	    delete infile;
+	    else if(ng_geometry!=nullptr) mesh->SetGeometry(ng_geometry);
 	  }),py::call_guard<py::gil_scoped_release>())
     // static_cast<void(Mesh::*)(const string & name)>(&Mesh::Load))
     .def("Save", static_cast<void(Mesh::*)(const string & name)const>(&Mesh::Save),py::call_guard<py::gil_scoped_release>())