Mesh reader extended to UGGrid, Filesystem-Library added, String utilities added

dune/amdis/CreatorInterface.hpp

@@ -3,6 +3,7 @@
 #include <memory>
 #include <string>
 
+#include <dune/amdis/common/Utility.hpp>
 #include <dune/amdis/Log.hpp>
 
 namespace AMDiS
@@ -30,7 +31,7 @@ namespace AMDiS
      * Must be implemented by sub classes of CreatorInterface.
      * Creates a new instance of the sub class of BaseClass.
      */
-    virtual std::shared_ptr<BaseClass> create() = 0; 
+    virtual shared_ptr<BaseClass> create() = 0; 
   };
 
   /**
@@ -44,9 +45,9 @@ namespace AMDiS
   {
   public:
     
-    virtual std::shared_ptr<BaseClass> create() override
+    virtual shared_ptr<BaseClass> create() override
     { 
-      AMDIS_ERROR_EXIT("This create method should not be called. Call create(string) instead!"); 
+      AMDIS_ERROR_EXIT("Should not be called. Call create(string) instead!"); 
       return {}; 
     };
     
@@ -55,7 +56,7 @@ namespace AMDiS
      * Creates a new instance of the sub class of BaseClass by passing a 
      * string to the constructor.
      */
-    virtual std::shared_ptr<BaseClass> create(std::string) = 0;
+    virtual shared_ptr<BaseClass> create(std::string) = 0;
   };
   
   

dune/amdis/FileWriter.hpp

@@ -45,8 +45,8 @@ namespace AMDiS
       dir = "output";
       Parameters::get(base + "->output directory", dir);
       
-      vtkWriter = std::make_shared<Dune::VTKWriter<MeshView>>(meshView);
-      vtkSeqWriter = std::make_shared<Dune::VTKSequenceWriter<MeshView>>(vtkWriter, filename, dir, "");
+      vtkWriter = make_shared<Dune::VTKWriter<MeshView>>(meshView);
+      vtkSeqWriter = make_shared<Dune::VTKSequenceWriter<MeshView>>(vtkWriter, filename, dir, "");
     }
 
     
@@ -125,8 +125,8 @@ namespace AMDiS
 
   private:
     MeshView const& meshView;
-    std::shared_ptr<Dune::VTKWriter<MeshView>> vtkWriter;
-    std::shared_ptr<Dune::VTKSequenceWriter<MeshView>> vtkSeqWriter;
+    shared_ptr<Dune::VTKWriter<MeshView>> vtkWriter;
+    shared_ptr<Dune::VTKSequenceWriter<MeshView>> vtkSeqWriter;
     
     std::vector<std::string> names;
     std::array<std::vector<double>, nComponents> data_vectors;

dune/amdis/Mesh.hpp

@@ -13,55 +13,68 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include "Initfile.hpp"
-#include "Log.hpp"
+#include <dune/grid/albertagrid/albertareader.hh>
+#include <dune/grid/io/file/gmshreader.hh>
+
+#include <dune/amdis/Initfile.hpp>
+#include <dune/amdis/Log.hpp>
+#include <dune/amdis/common/Filesystem.hpp>
 
 namespace AMDiS
 {
-  struct _albertagrid {};
-  struct _uggrid {};
-  struct _yaspgrid {};
-  struct _unknowngrid {};
-
-  template <class Grid>
-  struct MeshTagImpl
+  namespace tag 
   {
-    using type = _unknowngrid;
-  };
-  
-  template <class Grid>
-  using MeshTag = typename MeshTagImpl<Grid>::type;
+    struct albertagrid {};
+    struct uggrid {};
+    struct yaspgrid {};
+    struct unknowngrid {};
+    
+  } // end namespace tag
   
-  // specialization for some grid types from DUNE
-#if HAVE_ALBERTA
-  template <int dim, int dimworld>
-  struct MeshTagImpl<Dune::AlbertaGrid<dim, dimworld>>
+  namespace Impl
   {
-    using type = _albertagrid;
-  };
+    template <class Grid>
+    struct MeshTag
+    {
+      using type = tag::unknowngrid;
+    };
+    
+    // specialization for some grid types from DUNE
+#if HAVE_ALBERTA
+    template <int dim, int dimworld>
+    struct MeshTag<Dune::AlbertaGrid<dim, dimworld>>
+    {
+      using type = tag::albertagrid;
+    };
 #endif
-  
+    
 #if HAVE_UG
-  template <int dim>
-  struct MeshTagImpl<Dune::UGGrid<dim>>
-  {
-    using type = _uggrid;
-  };
+    template <int dim>
+    struct MeshTag<Dune::UGGrid<dim>>
+    {
+      using type = tag::uggrid;
+    };
 #endif
+    
+    template <int dim, class Coordinates>
+    struct MeshTag<Dune::YaspGrid<dim, Coordinates>>
+    {
+      using type = tag::yaspgrid;
+    };
+    
+  } // end namespace Impl
   
-  template <int dim, class Coordinates>
-  struct MeshTagImpl<Dune::YaspGrid<dim, Coordinates>>
-  {
-    using type = _yaspgrid;
-  };
+  template <class Grid>
+  using MeshTag_t = typename Impl::MeshTag<Grid>::type;
   
   
+  /// A creator class for meshes. Each mesh needs different way of initialization
   template <class Grid>
   class MeshCreator
   {
-    static std::unique_ptr<Grid> create(std::string meshName)
+    static unique_ptr<Grid> create(std::string meshName)
     {
-      AMDIS_ERROR_EXIT("Not yet implemented");
+      AMDIS_ERROR_EXIT("Creator not yet implemented for this mesh type.");
     }
   };
     
@@ -71,7 +84,7 @@ namespace AMDiS
   {
     using Grid = Dune::AlbertaGrid<dim, dimworld>;
     
-    static std::unique_ptr<Grid> create(std::string meshName)
+    static unique_ptr<Grid> create(std::string meshName)
     {
       std::string macro_filename = "";
       Parameters::get(meshName + "->macro file name", macro_filename);
@@ -79,18 +92,53 @@ namespace AMDiS
       // TODO: if filename_extension is ".2d" or ".3d" read it directly from file
       // otherwise use a factory method
       
-      return std::make_unique<Grid>(macro_filename);
+      return make_unique<Grid>(macro_filename);
     }
   };
 #endif
   
   
+#if HAVE_UG
+  template <int dim>
+  struct MeshCreator<Dune::UGGrid<dim>>
+  {
+    using Grid = Dune::UGGrid<dim>;
+    
+    static unique_ptr<Grid> create(std::string meshName)
+    {
+      
+      std::string filename = "";
+      Parameters::get(meshName + "->macro file name", filename);
+      
+      if (!filename.empty()) {
+        Path fn(filename);
+        auto ext = fn.extension();
+
+#if HAVE_ALBERTA
+        if (ext == "1d" || ext == "2d" || ext == "3d") {
+          Dune::GridFactory<Grid> factory;
+          Dune::AlbertaReader<Grid> reader; 
+          reader.readGrid(filename, factory);
+          return unique_ptr<Grid>{factory.createGrid()};
+        }
+#endif
+        if (ext == "msh") {
+          Dune::GmshReader<Grid> reader;
+          return unique_ptr<Grid>{reader.read(filename)};
+        }
+      } else {
+        AMDIS_ERROR_EXIT("Construction of UGGrid without filename not yet implemented!");
+      }
+    }
+  };
+#endif
+  
   template <int dim, class T>
   struct MeshCreator<Dune::YaspGrid<dim, Dune::EquidistantCoordinates<T,dim>>>
   {
     using Grid = Dune::YaspGrid<dim, Dune::EquidistantCoordinates<T,dim>>;
     
-    static std::unique_ptr<Grid> create(std::string meshName)
+    static unique_ptr<Grid> create(std::string meshName)
     {
       Dune::FieldVector<double, dim> L; L = 1.0;  // extension of the domain
       Parameters::get(meshName + "->dimension", L);
@@ -102,7 +150,7 @@ namespace AMDiS
       
       // TODO: add more parameters for yasp-grid (see constructor)
       
-      return std::make_unique<Grid>(L, s);
+      return make_unique<Grid>(L, s);
     }
   };
   
@@ -112,7 +160,7 @@ namespace AMDiS
   {
     using Grid = Dune::YaspGrid<dim, Dune::EquidistantOffsetCoordinates<T, dim>>;
     
-    static std::unique_ptr<Grid> create(std::string meshName)
+    static unique_ptr<Grid> create(std::string meshName)
     {
       Dune::FieldVector<double, dim> lowerleft;  // Lower left corner of the domain
       Dune::FieldVector<double, dim> upperright; // Upper right corner of the domain
@@ -124,7 +172,7 @@ namespace AMDiS
       
       // TODO: add more parameters for yasp-grid (see constructor)
       
-      return std::make_unique<Grid>(lowerleft, upperright, s);
+      return make_unique<Grid>(lowerleft, upperright, s);
     }
   };
 

dune/amdis/Operator.hpp

@@ -30,9 +30,9 @@ namespace AMDiS
     }
     
     template <class... Args>
-    static std::shared_ptr<Self> zot(Args&&... args)
+    static shared_ptr<Self> zot(Args&&... args)
     {
-      auto op = std::make_shared<Self>();
+      auto op = make_shared<Self>();
       op->addZOT(std::forward<Args>(args)...);
       return op;
     }
@@ -59,9 +59,9 @@ namespace AMDiS
     }
     
     template <class... Args>
-    static std::shared_ptr<Self> fot(Args&&... args)
+    static shared_ptr<Self> fot(Args&&... args)
     {
-      auto op = std::make_shared<Self>();
+      auto op = make_shared<Self>();
       op->addFOT(std::forward<Args>(args)...);
       return op;
     }
@@ -86,9 +86,9 @@ namespace AMDiS
     }
     
     template <class... Args>
-    static std::shared_ptr<Self> sot(Args&&... args)
+    static shared_ptr<Self> sot(Args&&... args)
     {
-      auto op = std::make_shared<Self>();
+      auto op = make_shared<Self>();
       op->addSOT(std::forward<Args>(args)...);
       return op;
     }
@@ -97,7 +97,7 @@ namespace AMDiS
     /// Calls \ref zot(), \ref for() or \ref sot(), depending on template 
     /// parameter \p Order.
     template <size_t Order, class... Args>
-    static std::shared_ptr<Self> create(Args&&... args)
+    static shared_ptr<Self> create(Args&&... args)
     {
       return create(index_<Order>{}, std::forward<Args>(args)...);
     }
@@ -180,19 +180,19 @@ namespace AMDiS
     Self& addSOTImpl(Term const& term, const index_<I>, const index_<J>);
     
     template <class... Args>
-    static std::shared_ptr<Self> create(index_<0>, Args&&... args)
+    static shared_ptr<Self> create(index_<0>, Args&&... args)
     {
       return zot(std::forward<Args>(args)...);
     }
     
     template <class... Args>
-    static std::shared_ptr<Self> create(index_<1>, Args&&... args)
+    static shared_ptr<Self> create(index_<1>, Args&&... args)
     {
       return fot(std::forward<Args>(args)...);
     }
     
     template <class... Args>
-    static std::shared_ptr<Self> create(index_<2>, Args&&... args)
+    static shared_ptr<Self> create(index_<2>, Args&&... args)
     {
       return sot(std::forward<Args>(args)...);
     }

dune/amdis/ProblemInstat.hpp

@@ -76,7 +76,7 @@ namespace AMDiS
     ProblemType& problemStat;
 
     /// Solution of the last timestep.
-    std::shared_ptr<SystemVectorType> oldSolution;
+    shared_ptr<SystemVectorType> oldSolution;
   };
 
   

dune/amdis/ProblemInstat.inc.hpp

@@ -43,7 +43,7 @@ namespace AMDiS
       
       // create oldSolution
       std::vector<std::string> componentNames(size, name + "_uOld");
-      oldSolution = std::make_shared<SystemVectorType>(*problemStat.getFeSpaces(), componentNames);
+      oldSolution = make_shared<SystemVectorType>(*problemStat.getFeSpaces(), componentNames);
     }
   }
 

dune/amdis/ProblemStat.hpp

@@ -102,12 +102,12 @@ namespace AMDiS
                            double* factor = NULL, 
                            double* estFactor = NULL);
 
-    void addMatrixOperator(std::shared_ptr<OperatorType> op,
+    void addMatrixOperator(shared_ptr<OperatorType> op,
                            int i, int j,
                            double* factor = NULL, 
                            double* estFactor = NULL);
     
-    void addMatrixOperator(std::map< std::pair<int,int>, std::shared_ptr<OperatorType> > ops);
+    void addMatrixOperator(std::map< std::pair<int,int>, shared_ptr<OperatorType> > ops);
     
     /// Adds an operator to \ref rhs.
     void addVectorOperator(OperatorType& op, 
@@ -115,12 +115,12 @@ namespace AMDiS
                            double* factor = NULL, 
                            double* estFactor = NULL);
     
-    void addVectorOperator(std::shared_ptr<OperatorType> op, 
+    void addVectorOperator(shared_ptr<OperatorType> op, 
                            int i,
                            double* factor = NULL, 
                            double* estFactor = NULL);
     
-    void addVectorOperator(std::map< int, std::shared_ptr<OperatorType> > ops);
+    void addVectorOperator(std::map< int, shared_ptr<OperatorType> > ops);
     
 
     /// Adds a Dirichlet boundary condition
@@ -219,7 +219,7 @@ namespace AMDiS
           "No mesh name specified for '" << name << "->mesh'!");
       
       mesh = MeshCreator<Mesh>::create(meshName);
-      meshView = std::make_shared<MeshView>(mesh->leafGridView());
+      meshView = make_shared<MeshView>(mesh->leafGridView());
       
       AMDIS_MSG("Create mesh:");
       AMDIS_MSG("#elements = "    << mesh->size(0));
@@ -230,16 +230,16 @@ namespace AMDiS
     
     void createFeSpaces()
     {
-      feSpaces = std::make_shared<FeSpaces>(constructTuple<FeSpaces>(*meshView));
+      feSpaces = make_shared<FeSpaces>(constructTuple<FeSpaces>(*meshView));
     }
     
     void createMatricesAndVectors()
     {
-      systemMatrix = std::make_shared<SystemMatrixType>(*feSpaces);
-      solution = std::make_shared<SystemVectorType>(*feSpaces, componentNames);
+      systemMatrix = make_shared<SystemMatrixType>(*feSpaces);
+      solution = make_shared<SystemVectorType>(*feSpaces, componentNames);
       
       auto rhsNames = std::vector<std::string>(nComponents, "rhs");
-      rhs = std::make_shared<SystemVectorType>(*feSpaces, rhsNames);
+      rhs = make_shared<SystemVectorType>(*feSpaces, rhsNames);
     }
     
     void createSolver()
@@ -255,7 +255,7 @@ namespace AMDiS
     
     void createFileWriter()
     {
-      filewriter = std::make_shared<FileWriter<Traits>>(name + "->output", 
+      filewriter = make_shared<FileWriter<Traits>>(name + "->output", 
                                                         *meshView, 
                                                         componentNames);
     }
@@ -271,13 +271,13 @@ namespace AMDiS
     bool getElementMatrix(RowView const& rowView,
 			  ColView const& colView,
 			  ElementMatrix& elementMatrix,
-			  std::list<std::shared_ptr<OperatorType>>& operators,
+			  std::list<shared_ptr<OperatorType>>& operators,
                           std::list<double*> const& factors);
     
     template <class RowView>
     bool getElementVector(RowView const& rowView,
 			  ElementVector& elementVector,
-			  std::list<std::shared_ptr<OperatorType>>& operators,
+			  std::list<shared_ptr<OperatorType>>& operators,
                           std::list<double*> const& factors);
     
     
@@ -311,49 +311,57 @@ namespace AMDiS
     std::vector<std::string> componentNames;
     
     /// Mesh of this problem.
-    std::shared_ptr<Mesh> mesh; // TODO: generalize to multi-mesh problems
+    shared_ptr<Mesh> mesh; // TODO: generalize to multi-mesh problems
     
     /// Name of the mesh
     std::string meshName;
     
     /// A gridView object
-    std::shared_ptr<MeshView> meshView;
+    shared_ptr<MeshView> meshView;
     
     /// Pointer to the meshes for the different problem components
     std::vector<Mesh*> componentMeshes;
     
     /// FE spaces of this problem.
-    std::shared_ptr<FeSpaces> feSpaces; // eventuell const
+    shared_ptr<FeSpaces> feSpaces; // eventuell const
     
     /// A FileWriter object
-    std::shared_ptr<FileWriter<Traits>> filewriter;
+    shared_ptr<FileWriter<Traits>> filewriter;
     
     /// An object of the linearSolver Interface
-    std::shared_ptr<LinearSolverType> linearSolver;
+    shared_ptr<LinearSolverType> linearSolver;
     
     /// A block-matrix that is filled during assembling
-    std::shared_ptr<SystemMatrixType> systemMatrix;
+    shared_ptr<SystemMatrixType> systemMatrix;
     
     /// A block-vector with the solution components
-    std::shared_ptr<SystemVectorType> solution;
+    shared_ptr<SystemVectorType> solution;
     
     /// A block-vector (load-vector) corresponding to the right.hand side 
     /// of the equation, filled during assembling
-    std::shared_ptr<SystemVectorType> rhs;
+    shared_ptr<SystemVectorType> rhs;
+    
+    
+    template <class T>
+    using MatrixEntries = std::map< std::pair<int,int>, std::list<T> >;
+    
+    template <class T>
+    using VectorEntries = std::map< int, std::list<T> >;
+  
     
     /// A map (i,j) -> list<DirichleBC> string a boundary conditions for 
     /// each matrix block
-    IdxPairList<DirichletBC<WorldVector>> dirichletBc;
+    MatrixEntries<shared_ptr<DirichletBC<WorldVector>>> dirichletBc;
     
     /// A map (i,j) -> list<OperatorType> string the differential operators for 
     /// each matrix block
-    IdxPairList<OperatorType> matrixOperators;
-    std::map<std::pair<int,int>, std::list<double*>> matrixFactors;
+    MatrixEntries<shared_ptr<OperatorType>> matrixOperators;
+    MatrixEntries<double*>                  matrixFactors;
     
     /// A map (i) -> list<OperatorType> string the differential operators for 
     /// each rhs block
-    IdxList<OperatorType> vectorOperators;
-    std::map<int, std::list<double*>> vectorFactors;
+    VectorEntries<shared_ptr<OperatorType>> vectorOperators;
+    VectorEntries<double*>                  vectorFactors;
   };
   
   

dune/amdis/ProblemStat.inc.hpp

@@ -107,7 +107,7 @@ namespace AMDiS
 
 
   template <class Traits>
-  void ProblemStatSeq<Traits>::addMatrixOperator(std::shared_ptr<OperatorType> op, 
+  void ProblemStatSeq<Traits>::addMatrixOperator(shared_ptr<OperatorType> op, 
                                                  int i, int j,
                                                  double* factor, 
                                                  double* estFactor)
@@ -126,18 +126,18 @@ namespace AMDiS
                                                  double* factor, 
                                                  double* estFactor)
   {
-    addMatrixOperator(std::make_shared<OperatorType>(op), i, j, factor, estFactor);
+    addMatrixOperator(make_shared<OperatorType>(op), i, j, factor, estFactor);
   }
 
   template <class Traits>
-  void ProblemStatSeq<Traits>::addMatrixOperator(std::map< std::pair<int,int>, std::shared_ptr<OperatorType> > ops)
+  void ProblemStatSeq<Traits>::addMatrixOperator(std::map< std::pair<int,int>, shared_ptr<OperatorType> > ops)
   {
     for (auto op : ops)
       addMatrixOperator(op.second, op.first.first, op.first.second);
   }
   
   template <class Traits>
-  void ProblemStatSeq<Traits>::addVectorOperator(std::shared_ptr<OperatorType> op, 
+  void ProblemStatSeq<Traits>::addVectorOperator(shared_ptr<OperatorType> op, 
                                                  int i,
                                                  double* factor, 
                                                  double* estFactor)
@@ -155,12 +155,12 @@ namespace AMDiS
                                                  double* factor, 
                                                  double* estFactor)
   {
-    addVectorOperator(std::make_shared<OperatorType>(op), i, factor, estFactor);
+    addVectorOperator(make_shared<OperatorType>(op), i, factor, estFactor);
   }
 
   
   template <class Traits>
-  void ProblemStatSeq<Traits>::addVectorOperator(std::map< int, std::shared_ptr<OperatorType> > ops)
+  void ProblemStatSeq<Traits>::addVectorOperator(std::map< int, shared_ptr<OperatorType> > ops)
   {
     for (auto op : ops)
       addVectorOperator(op.second, op.first);
@@ -360,7 +360,7 @@ namespace AMDiS
   bool ProblemStatSeq<Traits>::getElementMatrix(RowView const& rowLocalView,
                                                 ColView const& colLocalView,
                                                 ElementMatrix& elementMatrix,
-                                                std::list<std::shared_ptr<OperatorType>>& operators,
+                                                std::list<shared_ptr<OperatorType>>& operators,
                                                 std::list<double*> const& factors)
   {
     const auto nRows = rowLocalView.tree().finiteElement().size();
@@ -390,7 +390,7 @@ namespace AMDiS
     template <class RowView>
   bool ProblemStatSeq<Traits>::getElementVector(RowView const& rowLocalView,
                                                 ElementVector& elementVector,
-                                                std::list<std::shared_ptr<OperatorType>>& operators,
+                                                std::list<shared_ptr<OperatorType>>& operators,
                                                 std::list<double*> const& factors)
   {
     const auto nRows = rowLocalView.tree().finiteElement().size();

dune/amdis/common/Filesystem.hpp

+#pragma once
+
+#include <string>
+#include <stdio.h>  /* defines FILENAME_MAX */
+#include <cctype>
+
+#ifdef WINDOWS
+    #include <direct.h>
+    #define getCurrentDir _getcwd
+#else
+    #include <unistd.h>
+    #define getCurrentDir getcwd
+#endif
+ 
+#include <dune/amdis/common/String.hpp> 
+
+namespace AMDiS
+{
+  // A minimalistic filesystem class
+  class Path
+  {
+  public:
+    Path(std::string p)
+    {
+      split(p, pathParts);
+      path_ = normalize(pathParts);
+    }
+    
+    // return filename part of path
+    std::string filename() const
+    {
+      return pathParts.back();
+    }
+    
+    // get filename extension
+    std::string extension() const
+    {
+      auto fn = filename();
+      auto pos = fn.find_last_of('.');
+      return fn.substr(pos+1);
+    }
+    
+    // get current working directory
+    std::string currentPath() const
+    {
+      char cwd_[FILENAME_MAX];
+      getCurrentDir(cwd_, sizeof(cwd_));
+      std::string cwd(cwd_);
+      return trim(cwd);
+    }