diff --git a/src/amdis/DataTransfer.inc.hpp b/src/amdis/DataTransfer.inc.hpp
index 2e52110d3fa96b1759cedceaa6cd3f411fc9281e..64319984308b8a81e43fab67ddcba599a4be8bcb 100644
--- a/src/amdis/DataTransfer.inc.hpp
+++ b/src/amdis/DataTransfer.inc.hpp
@@ -365,7 +365,7 @@ restrictLocal(Element const& father, NodeElementData& fatherDOFs, Trafo const& t
         for (std::size_t i = 0; i < shapeValues.size(); ++i)
           y += shapeValues[i] * childDOFs[i];
 
-        fatherDOFsTemp_[currentDOF] = y;
+        fatherDOFsTemp_[currentDOF] = T(y);
         finishedDOFs_[currentDOF++] = true;
         return y;
       }
diff --git a/src/amdis/common/CMakeLists.txt b/src/amdis/common/CMakeLists.txt
index 006a5bc0db252077cf651779bedf8a8f0b7e6c0d..f14c1c2468855aab773327cc6784d7b3a23de38c 100644
--- a/src/amdis/common/CMakeLists.txt
+++ b/src/amdis/common/CMakeLists.txt
@@ -27,6 +27,7 @@ install(FILES
     MultiTypeVector.hpp
     Range.hpp
     Resize.hpp
+    QuadMath.hpp
     StaticSize.hpp
     String.hpp
     SwitchCases.hpp
diff --git a/src/amdis/common/FieldMatVec.inc.hpp b/src/amdis/common/FieldMatVec.inc.hpp
index dc4fe7fbcdc167889c082725a8f6fba6e23b4f15..dc68bebff5019a71e5f2bc819d76a572ec906ed3 100644
--- a/src/amdis/common/FieldMatVec.inc.hpp
+++ b/src/amdis/common/FieldMatVec.inc.hpp
@@ -346,13 +346,15 @@ auto two_norm(T const& x)
 template <class T, int N>
 auto two_norm(FieldVector<T, N> const& x)
 {
-  return std::sqrt(unary_dot(x));
+  using std::sqrt;
+  return sqrt(unary_dot(x));
 }
 
 template <class T, int N>
 auto two_norm(FieldMatrix<T, 1, N> const& x)
 {
-  return std::sqrt(unary_dot(x));
+  using std::sqrt;
+  return sqrt(unary_dot(x));
 }
 
 /** \ingroup vector_norms
diff --git a/src/amdis/common/QuadMath.hpp b/src/amdis/common/QuadMath.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..1c0ba7e08203935a0047233724a5cc4f810ed1ee
--- /dev/null
+++ b/src/amdis/common/QuadMath.hpp
@@ -0,0 +1,115 @@
+#pragma once
+
+#if HAVE_QUADMATH
+#include <functional>
+#include <type_traits>
+
+#include <dune/common/hash.hh>
+#include <dune/common/quadmath.hh>
+
+#include <amdis/common/DerivativeTraits.hpp>
+#include <amdis/common/StaticSize.hpp>
+#include <amdis/common/ValueCategory.hpp>
+#include <amdis/gridfunctions/ConstantGridFunction.hpp>
+
+namespace std
+{
+  template <class T>
+  struct common_type<Dune::Float128, T>
+  {
+    using type = Dune::Float128;
+  };
+
+  template <class T>
+  struct common_type<T, Dune::Float128>
+  {
+    using type = Dune::Float128;
+  };
+
+  template <>
+  struct common_type<Dune::Float128, Dune::Float128>
+  {
+    using type = Dune::Float128;
+  };
+
+  template <>
+  struct hash<Dune::Float128>
+  {
+    typedef Dune::Float128 argument_type;
+    typedef std::size_t result_type;
+
+    std::size_t operator()(const Dune::Float128& arg) const
+    {
+      hash<long double> hasher_ld;
+      return hasher_ld((long double)(arg));
+    }
+  };
+
+  template <>
+  struct hash<const Dune::Float128>
+      : public hash<Dune::Float128>
+  {};
+
+} // end namespace std
+
+
+namespace Dune
+{
+  namespace Impl
+  {
+    // specialization for float arguments due to ambiguity
+    template <class T,
+      std::enable_if_t<not std::is_integral<T>::value && std::is_arithmetic<T>::value, int> = 0>
+    inline Float128 pow(const Float128& x, const T& p)
+    {
+      return powq(float128_t(x), float128_t(p));
+    }
+
+  } // end namespace Impl
+} // end namespace Dune
+
+
+namespace AMDiS
+{
+  namespace Concepts
+  {
+    namespace Definition
+    {
+      template <>
+      struct ConstantToGridFunction<Dune::Float128>
+        : std::true_type {};
+
+    } // end namespace Definition
+  } // end namespace Concepts
+
+  namespace Impl
+  {
+    template <>
+    struct SizeImpl<Dune::Float128>
+        : std::integral_constant<std::size_t, 1> {};
+
+    template <>
+    struct RowsImpl<Dune::Float128>
+        : std::integral_constant<std::size_t, 1> {};
+
+    template <>
+    struct ColsImpl<Dune::Float128>
+        : std::integral_constant<std::size_t, 1> {};
+
+  } // end namespace Impl
+
+  template <class K, int N>
+  struct DerivativeTraits<Dune::Float128(Dune::FieldVector<K,N>), tag::gradient>
+  {
+    using Range = Dune::FieldVector<Dune::Float128,N>;
+  };
+
+  template <>
+  struct ValueCategory<Dune::Float128>
+  {
+    using type = tag::scalar;
+  };
+
+} // end namespace AMDiS
+
+#endif // HAVE_QUADMATH
diff --git a/src/amdis/functions/GlobalIdSet.hpp b/src/amdis/functions/GlobalIdSet.hpp
index 29f56a69292cacd7817973faf066f38615a6134a..ed3e7a8b62f3ba3d169ed046e19c997e4124c224 100644
--- a/src/amdis/functions/GlobalIdSet.hpp
+++ b/src/amdis/functions/GlobalIdSet.hpp
@@ -29,9 +29,6 @@ namespace Dune
     template <class MI, class IMS, class... SPB>
     class CompositePreBasis;
 
-    template <class GV, int k, class MI>
-    class LagrangePreBasis;
-
     template <class MI, class IMS, class SPB, std::size_t C>
     class PowerPreBasis;
 
diff --git a/src/amdis/linearalgebra/petsc/Traits.hpp b/src/amdis/linearalgebra/petsc/Traits.hpp
index bba67a939d2cc5c79825ffcefc7605b2439f7e19..fb5e1598eab595a38dc032b861e94d5678c639de 100644
--- a/src/amdis/linearalgebra/petsc/Traits.hpp
+++ b/src/amdis/linearalgebra/petsc/Traits.hpp
@@ -16,7 +16,7 @@ namespace AMDiS
   template <class Basis, class T = double>
   struct BackendTraits
   {
-    static_assert(std::is_convertible<T,PetscScalar>::value, "");
+    static_assert(std::is_same<T,PetscScalar>::value, "");
 
     using Mat = ::Mat;
     using Vec = ::Vec;
diff --git a/test/ProblemStatTest.cpp b/test/ProblemStatTest.cpp
index 8c6b36dbbadcb315ce8186ae9307eb7a06df5b05..b44b52f46ec443cece066c846c2bf38248c7c56a 100644
--- a/test/ProblemStatTest.cpp
+++ b/test/ProblemStatTest.cpp
@@ -2,20 +2,33 @@
 #include <amdis/AdaptStationary.hpp>
 #include <amdis/LocalOperators.hpp>
 #include <amdis/ProblemStat.hpp>
+#include <amdis/common/QuadMath.hpp>
 
-using namespace AMDiS;
-
-using Grid = Dune::YaspGrid<2>;
+#include <dune/common/version.hh>
+#include <dune/grid/yaspgrid.hh>
 
-template <class T>
-struct Param
-    : DefaultBasisCreator<Grid, Impl::LagrangePreBasisCreator<1,1>, T>
-{};
+using namespace AMDiS;
 
 template <class T>
 void test()
 {
-  ProblemStat<Param<T>> prob("ellipt");
+  // use T as coordinate type
+  using Grid = Dune::YaspGrid<2, Dune::EquidistantCoordinates<T,2>>;
+  Grid grid({T(1), T(1)}, {8,8});
+
+  // use T as range type for basis
+  using namespace Dune::Functions::BasisFactory;
+#if DUNE_VERSION_LT(DUNE_FUNCTIONS,2,7)
+  auto basis = makeBasis(grid.leafGridView(), power<1>(lagrange<1>(), flatLexicographic()));
+#else
+  auto basis = makeBasis(grid.leafGridView(), power<1>(lagrange<1,T>(), flatLexicographic()));
+#endif
+  using Basis = decltype(basis);
+
+  // use T as coefficient type
+  using Param = DefaultProblemTraits<Basis, T>;
+
+  ProblemStat<Param> prob("ellipt", grid, basis);
   prob.initialize(INIT_ALL);
   prob.boundaryManager()->setBoxBoundary({1,1,2,2});
 
@@ -45,5 +58,9 @@ int main(int argc, char** argv)
   test<double>();
   test<long double>();
 
+#if HAVE_QUADMATH
+  test<Dune::Float128>();
+#endif
+
   return 0;
 }