From 6f6d1bcb4d4f76de32b97788c9374d5ab9060f34 Mon Sep 17 00:00:00 2001
From: Andreas Dedner <a.s.dedner@warwick.ac.uk>
Date: Wed, 3 Oct 2012 22:24:50 +0100
Subject: [PATCH] A slightly simpler version for seriel and parallel call to
 refinement and some more debugging

---
 src/duneinterface/gitter_dune_pll_impl.cc | 166 +++++++++++++---------
 src/parallel/gitter_pll_sti.cc            |   4 +-
 src/serial/gitter_sti.cc                  |  77 +++++-----
 src/serial/gitter_sti.h                   |   2 -
 4 files changed, 145 insertions(+), 104 deletions(-)

diff --git a/src/duneinterface/gitter_dune_pll_impl.cc b/src/duneinterface/gitter_dune_pll_impl.cc
index 4317ba58f..161df1a48 100644
--- a/src/duneinterface/gitter_dune_pll_impl.cc
+++ b/src/duneinterface/gitter_dune_pll_impl.cc
@@ -1416,27 +1416,40 @@ void GitterDunePll :: tovtk( const std::string &fn )
   VertexList vertexList;
 
   int nCells = 0;
-
+  int nBnd = 0;
+#if 0
+  typedef LeafIterator < Gitter::helement_STI > Iterator;
+  Iterator w (*this) ;
+  typedef LeafIterator < Gitter::hbndseg_STI > BndIterator;
+  BndIterator wbnd (*this) ;
+#else
+  typedef LevelIterator < Gitter::helement_STI > Iterator;
+  Iterator w (*this,0) ;
+  typedef LevelIterator < Gitter::hbndseg_STI > BndIterator;
+  // typedef LevelIterator < Gitter::hface_STI > BndIterator;
+  BndIterator wbnd (*this,0) ;
+#endif
   // loop to find vertexList and count cells
   {
     typedef Objects :: tetra_IMPL tetra_IMPL ;
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
-      {
-	tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
-	for (int i=0;i<4;++i)
-	  {
-	    Vertex v(3);
-	    for (int k=0;k<3;++k) 
-	      v[k] = item->myvertex(i)->Point()[k];
-	    vertexList[ item->myvertex(i)->getIndex() ] = make_pair(-1,v);
-	  }
-	++nCells;
-      }
+    {
+	    tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
+	    for (int i=0;i<4;++i)
+	    {
+	      Vertex v(3);
+	      for (int k=0;k<3;++k) 
+	        v[k] = item->myvertex(i)->Point()[k];
+	      vertexList[ item->myvertex(i)->getIndex() ] = make_pair(-1,v);
+	    }
+	    ++nCells;
+    }
+    for (wbnd->first () ; ! wbnd->done () ; wbnd->next ())
+      ++nBnd;
   }
 
   vtuFile << "    <Piece NumberOfPoints=\"" << vertexList.size() << "\" "
-	  << "NumberOfCells=\"" << nCells << "\">" << std::endl;
+	        << "NumberOfCells=\"" << nCells+nBnd << "\">" << std::endl;
 
   // cell data
   {
@@ -1445,12 +1458,22 @@ void GitterDunePll :: tovtk( const std::string &fn )
     vtuFile << "          ";
 
     typedef Objects :: tetra_IMPL tetra_IMPL ;
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
-      {
-	tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
-	vtuFile << item->getIndex() << " ";
-      }
+    {
+	    tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
+	    // vtuFile << item->getIndex() << " ";
+      bool ok = true;
+      for (int k=0;k<4;++k)
+        ok &= item->myneighbour( k ).first->isRealObject();
+      vtuFile << ((ok)?1:2) << " ";
+    }
+
+    vtuFile << std::endl;
+    for (wbnd->first () ; ! wbnd->done () ; wbnd->next ())
+    {
+      vtuFile << -1 << " ";
+    }
+
 
     vtuFile << std::endl;
     vtuFile << "        </DataArray>" << std::endl;
@@ -1482,15 +1505,23 @@ void GitterDunePll :: tovtk( const std::string &fn )
     vtuFile << "         ";
 
     typedef Objects :: tetra_IMPL tetra_IMPL ;
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
-      {
-	tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
-	for (int i=0;i<4;++i)
-	  {
-	    vtuFile << " " << vertexList[item->myvertex(i)->getIndex()].first;
-	  }
-      }
+    {
+	    tetra_IMPL* item = ((tetra_IMPL *) &w->item ());
+	    for (int i=0;i<4;++i)
+	    {
+	      vtuFile << " " << vertexList[item->myvertex(i)->getIndex()].first;
+	    }
+    }
+    for (wbnd->first () ; ! wbnd->done () ; wbnd->next ())
+    {
+	    hbndseg3_GEO* item = ((hbndseg3_GEO *) &wbnd->item ());
+	    // hface3_GEO* item = ((hface3_GEO *) &wbnd->item ());
+	    for (int i=0;i<3;++i)
+	    {
+	      vtuFile << " " << vertexList[item->myvertex(0,i)->getIndex()].first;
+	    }
+    }
     vtuFile << std::endl;
     vtuFile << "        </DataArray>" << std::endl;
 
@@ -1499,9 +1530,13 @@ void GitterDunePll :: tovtk( const std::string &fn )
     vtuFile << "         ";
 
     for( int i = 0; i < nCells; ++i )
-      {
-	vtuFile << " " << (i+1)*4;
-      }
+    {
+	    vtuFile << " " << (i+1)*4;
+    }
+    for( int i = 0; i < nBnd; ++i )
+    {
+	    vtuFile << " " << nCells*4 + (i+1)*3;
+    }
     vtuFile << std::endl;
 
     vtuFile << "        </DataArray>" << std::endl;
@@ -1511,9 +1546,13 @@ void GitterDunePll :: tovtk( const std::string &fn )
     vtuFile << "         ";
 
     for( int i = 0; i < nCells; ++i )
-      {
-	vtuFile << " " << 10; // 10 for tetrahedra
-      }
+    {
+	    vtuFile << " " << 10; // 10 for tetrahedra
+    }
+    for( int i = 0; i < nBnd; ++i )
+    {
+	    vtuFile << " " << 5; // 5 for triangle
+    }
     vtuFile << std::endl;
 
     vtuFile << "        </DataArray>" << std::endl;
@@ -1527,38 +1566,37 @@ void GitterDunePll :: tovtk( const std::string &fn )
   vtuFile.close();
 
   if( myrank == 0 )
-    {
-      std::ostringstream pllss;
-
-      if( fn.substr(fn.find_last_of(".") + 1) == "vtu" )
-	{
-	  pllss << fn.substr(0,fn.find_last_of(".") + 1) << "pvtu";
-	}
-      else
-	{
-	  pllss << fn << ".pvtu";
-	}
-      std::ofstream pvtuFile;
-      pvtuFile.open( pllss.str().c_str() );
+  {
+    std::ostringstream pllss;
+    if( fn.substr(fn.find_last_of(".") + 1) == "vtu" )
+	  {
+	    pllss << fn.substr(0,fn.find_last_of(".") + 1) << "pvtu";
+	  }
+    else
+	  {
+	    pllss << fn << ".pvtu";
+	  }
+    std::ofstream pvtuFile;
+    pvtuFile.open( pllss.str().c_str() );
   
-      pvtuFile << "<?xml version=\"1.0\"?>" << std::endl;
-      pvtuFile << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
-      pvtuFile << "  <PUnstructuredGrid GhostLevel=\"0\">" << std::endl;
-      pvtuFile << "    <PCellData Scalars=\"cell-id\">" << std::endl;
-      pvtuFile << "      <PDataArray type=\"Float32\" Name=\"cell-id\" />" << std::endl;
-      pvtuFile << "    </PCellData>" << std::endl;
-      pvtuFile << "    <PPoints>" << std::endl;
-      pvtuFile << "      <PDataArray type=\"Float32\" NumberOfComponents=\"3\" />" << std::endl;
-      pvtuFile << "    </PPoints>" << std::endl;
-      for( int p = 0; p < nProc; ++p )
-	pvtuFile << "    <Piece Source=\"p" << p << "-" << fn << "\" />" << std::endl;
-      pvtuFile << "  </PUnstructuredGrid>" << std::endl;
-      pvtuFile << "</VTKFile>" << std::endl;
-
-      pvtuFile.close();
-
-      std::cout << "data written to " << pllss.str() << std::endl;
-    }
+    pvtuFile << "<?xml version=\"1.0\"?>" << std::endl;
+    pvtuFile << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
+    pvtuFile << "  <PUnstructuredGrid GhostLevel=\"0\">" << std::endl;
+    pvtuFile << "    <PCellData Scalars=\"cell-id\">" << std::endl;
+    pvtuFile << "      <PDataArray type=\"Float32\" Name=\"cell-id\" />" << std::endl;
+    pvtuFile << "    </PCellData>" << std::endl;
+    pvtuFile << "    <PPoints>" << std::endl;
+    pvtuFile << "      <PDataArray type=\"Float32\" NumberOfComponents=\"3\" />" << std::endl;
+    pvtuFile << "    </PPoints>" << std::endl;
+    for( int p = 0; p < nProc; ++p )
+	    pvtuFile << "    <Piece Source=\"p" << p << "-" << fn << "\" />" << std::endl;
+    pvtuFile << "  </PUnstructuredGrid>" << std::endl;
+    pvtuFile << "</VTKFile>" << std::endl;
+
+    pvtuFile.close();
+
+    std::cout << "data written to " << pllss.str() << std::endl;
+  }
 }
 
 
diff --git a/src/parallel/gitter_pll_sti.cc b/src/parallel/gitter_pll_sti.cc
index 225f991bd..4ed705c86 100644
--- a/src/parallel/gitter_pll_sti.cc
+++ b/src/parallel/gitter_pll_sti.cc
@@ -868,7 +868,7 @@ bool GitterPll :: adapt () {
            << _refineLoops << ") " << u2 << " " << u3 << endl ;
     }
     notifyGridChanges () ;
-    // loadBalancerGridChangesNotify () ;
+    loadBalancerGridChangesNotify () ;
     conformClosure = mpAccess().gmax( !markNonConform() );
   } while (conformClosure);
   ++adaptstep;
@@ -1167,7 +1167,9 @@ void GitterPll :: loadBalancerGridChangesNotify ()
 
     if( ldbMth ) 
     {
+      tovtk("pre.vtu");
       repartitionMacroGrid (db) ;
+      tovtk("post.vtu");
       notifyMacroGridChanges () ;
     }
   }
diff --git a/src/serial/gitter_sti.cc b/src/serial/gitter_sti.cc
index 85346e92f..a9d4c5200 100644
--- a/src/serial/gitter_sti.cc
+++ b/src/serial/gitter_sti.cc
@@ -293,49 +293,41 @@ void Gitter :: printsize () {
 
 int nr = 0;
 int adaptstep = 0;
-bool Gitter :: doRefine ( const bool conformingClosure ) 
+bool Gitter :: refine () 
 {
   assert (debugOption (20) ? (cout << "**INFO GitterDuneBasis :: refine ()" << endl, 1) : 1) ;
   bool x = true ;
   leaf_element__macro_element__iterator i (container ()) ;
-  do
-  {
-    x = true ;
-    // refine marked elements
-    for( i.first(); ! i.done() ; i.next()) x &= i.item ().refine () ;
-
-    if( conformingClosure ) 
-    {
-      // check for conformity
-      // if noconform break;
-      std::cout << "check non conform refinement" << std::endl;
-      x = true ;
-      for( i.first(); ! i.done() ; i.next()) { x &= i.item ().markNonConform () ; }
-      std::ostringstream ss;
-      int filenr = adaptstep*100+nr;
-      ss << "ref-" << ZeroPadNumber(filenr) << ".vtu";
-      tovtk(  ss.str() );
-      ++nr;
-    }
-    // break;
-    if (x) break;
-  }
-  while (1);  // need something here on required conformity
-  ++adaptstep;
+  // refine marked elements
+  for( i.first(); ! i.done() ; i.next()) x &= i.item ().refine () ;
+  std::ostringstream ss;
+  int filenr = adaptstep*100+nr;
+  ss << "ref-" << ZeroPadNumber(filenr) << ".vtu";
+  tovtk(  ss.str() );
+  ++nr;
   return  x;
 }
 
-// refine without conforming closure 
+/*
+int nr = 0;
+int adaptstep = 0;
 bool Gitter :: refine () 
 {
-  return doRefine( false );
-}
-
-// refine including iterations for conforming closure 
-bool Gitter :: refineConforming ()
-{
-  return doRefine( true );
+  assert (debugOption (20) ? (cout << "**INFO GitterDuneBasis :: refine ()" << endl, 1) : 1) ;
+  bool x = true ;
+  leaf_element__macro_element__iterator i (container ()) ;
+  for( i.first(); ! i.done() ; i.next() ) 
+  {
+    x &= i.item ().refine () ;
+  }
+	std::ostringstream ss;
+  int filenr = adaptstep*1000+nr;
+	ss << "ref-" << ZeroPadNumber(filenr) << ".vtu";
+  tovtk(  ss.str() );
+  ++nr;
+  return  x;
 }
+*/
 
 bool Gitter :: markNonConform()
 {
@@ -377,10 +369,12 @@ void Gitter :: tovtk( const std::string &fn )
   int nCells = 0;
 
   typedef Gitter :: Geometric :: tetra_GEO tetra_GEO ;
+  // typedef LeafIterator < Gitter::helement_STI > Iterator;
+  typedef LevelIterator < Gitter::helement_STI > Iterator;
+  Iterator w (*this,0) ;
 
   // loop to find vertexList and count cells
   {
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
     {
       tetra_GEO* item = ((tetra_GEO *) &w->item ());
@@ -404,7 +398,6 @@ void Gitter :: tovtk( const std::string &fn )
     vtuFile << "        <DataArray type=\"Float32\" Name=\"cell-id\" NumberOfComponents=\"1\">" << std::endl;
     vtuFile << "          ";
 
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
     {
       vtuFile << w->item ().getIndex() << " ";
@@ -439,7 +432,6 @@ void Gitter :: tovtk( const std::string &fn )
     vtuFile << "        <DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\">" << std::endl;
     vtuFile << "         ";
 
-    LeafIterator < Gitter::helement_STI > w (*this) ;
     for (w->first () ; ! w->done () ; w->next ())
     {
       tetra_GEO* item = ((tetra_GEO *) &w->item ());
@@ -491,7 +483,18 @@ bool Gitter :: adapt ()
   assert (! iterators_attached ()) ;
   const int start = clock () ;
 
-  bool refined = refineConforming ();
+  bool x;
+  bool refined = true;
+  do {
+    refined &= refine ();
+    // check for conformity
+    // if noconform break;
+    std::cout << "check non conform refinement" << std::endl;
+    x = markNonConform();
+  }
+  while (!x);  // need something here on required conformity
+  ++adaptstep;
+
   if (!refined) {
     cerr << "**WARNUNG (IGNORIERT) Verfeinerung nicht vollst\"andig (warum auch immer)\n" ;
     cerr << "  diese Option ist eigentlich dem parallelen Verfeinerer vorbehalten.\n" ;
diff --git a/src/serial/gitter_sti.h b/src/serial/gitter_sti.h
index 22bc1a522..0223aae26 100644
--- a/src/serial/gitter_sti.h
+++ b/src/serial/gitter_sti.h
@@ -2021,8 +2021,6 @@ public:
 protected :
   // methods for refining and coarsening
   virtual bool refine () ;
-  virtual bool refineConforming () ;
-  bool doRefine( const bool conformingClosure ) ;
   virtual bool markNonConform () ;
   virtual void coarse () ;
 
-- 
GitLab