diff --git a/CHANGELOG.md b/CHANGELOG.md
index 5d44ec767..9f5210a00 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -52,6 +52,8 @@ SPDX-License-Identifier: LGPL-3.0-or-later
   ([#304](https://github.com/ikarus-project/ikarus/pull/304))
     - This can be used, for instance, to apply concentrated forces or to add spring stiffness in a particular direction.
     - Furthermore, a helper function to get the global index of a Lagrange node at the given global position is added.
+- Rework Python Interface for `DirichletValues` plus adding support to easily fix boundary DOFs of subspacebasis in C++ and Python ([#305](https://github.com/ikarus-project/ikarus/pull/305))
+- Rework the Python Interface for `DirichletValues` plus add support to easily fix boundary DOFs of `Subspacebasis` in C++ and Python ([#305](https://github.com/ikarus-project/ikarus/pull/305))
 
 ## Release v0.4 (Ganymede)
 
diff --git a/docs/website/01_framework/dirichletBCs.md b/docs/website/01_framework/dirichletBCs.md
index 0db59b14f..22cd30457 100644
--- a/docs/website/01_framework/dirichletBCs.md
+++ b/docs/website/01_framework/dirichletBCs.md
@@ -24,23 +24,25 @@ The interface of `#!cpp Ikarus::DirichletValues` is represented by the following
 Ikarus::DirichletValues dirichletValues2(basis); // (1)!
 void fixBoundaryDOFs(f); // (2)!
 void fixDOFs(f); // (3)!
-void fixIthDOF(i); // (4)!
+void setSingleDOF(i, flag); // (4)!
 const auto& basis() const; // (5)!
-bool isConstrained(std::size_t i) const; // (6)!
+bool isConstrained(i) const; // (6)!
 auto fixedDOFsize() const; // (7)!
 auto size() const ; // (8)!
+auto reset(); // (9)!
 ```
 
 1. Create class by inserting a global basis, [@sander2020dune] Chapter 10.
-2. Accepts a functor to fix boundary degrees of freedom. `f` is  a functor that will be called with the boolean vector of fixed boundary.
+2. Accepts a functor to fix boundary degrees of freedom. `f` is  a functor that will be called with the Boolean vector of fixed boundary.
  degrees of freedom and the usual arguments of `Dune::Functions::forEachBoundaryDOF`, as defined on page 388 of the Dune
    [@sander2020dune] book.
-3. A more general version of `fixBoundaryDOFs`. Here, a functor is to be provided that accepts a basis and the corresponding boolean
-4. A function that helps to fix the $i$-th degree of freedom
+3. A more general version of `fixBoundaryDOFs`. Here, a functor is to be provided that accepts a basis and the corresponding Boolean
+4. A function that helps to fix or unfix the $i$-th degree of freedom
    vector considering the Dirichlet degrees of freedom.
 5. Returns the underlying basis.
-6. Indicates whether the degree of freedom `i` is fixed.
+6. Indicates whether the degree of freedom $i$ is fixed.
 7. Returns the number of fixed degrees of freedom.
-8. Returns the number of all dirichlet degrees of freedom.
+8. Returns the number of all Dirichlet degrees of freedom.
+9. Resets the whole container
 
 \bibliography
diff --git a/ikarus/python/dirichletvalues/dirichletvalues.hh b/ikarus/python/dirichletvalues/dirichletvalues.hh
index 2b6320457..28ff9cf7b 100644
--- a/ikarus/python/dirichletvalues/dirichletvalues.hh
+++ b/ikarus/python/dirichletvalues/dirichletvalues.hh
@@ -8,11 +8,16 @@
 
 #pragma once
 
+#include <cstdlib>
+#include <string>
+
+#include "dune/common/classname.hh"
 #include <dune/functions/functionspacebases/lagrangebasis.hh>
 #include <dune/functions/functionspacebases/powerbasis.hh>
 #include <dune/grid/yaspgrid.hh>
 #include <dune/python/common/typeregistry.hh>
 #include <dune/python/functions/globalbasis.hh>
+#include <dune/python/functions/subspacebasis.hh>
 #include <dune/python/pybind11/eigen.h>
 #include <dune/python/pybind11/functional.h>
 #include <dune/python/pybind11/pybind11.h>
@@ -20,22 +25,116 @@
 #include <dune/python/pybind11/stl_bind.h>
 
 #include <ikarus/finiteelements/ferequirements.hh>
+
 // PYBIND11_MAKE_OPAQUE(std::vector<bool>);
 namespace Ikarus::Python {
 
+namespace Impl {
+  using FixBoundaryDOFsWithGlobalIndexFunction = std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int)>;
+
+  template <typename LV>
+  using FixBoundaryDOFsWithLocalViewFunction = std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int, LV&)>;
+
+  template <typename LV, typename IS>
+  using FixBoundaryDOFsWithIntersectionFunction =
+      std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int, LV&, const IS&)>;
+
+  template <typename Basis>
+  auto registerSubSpaceLocalView() {
+    pybind11::module scopedf = pybind11::module::import("dune.functions");
+    using LocalViewWrapper   = Dune::Python::LocalViewWrapper<Basis>;
+
+    auto includes = Dune::Python::IncludeFiles{"dune/python/functions/globalbasis.hh"};
+
+    Dune::Python::insertClass<Basis>(scopedf, "SubspaceBasis_" + Dune::className<typename Basis::PrefixPath>(),
+                                     Dune::Python::GenerateTypeName(Dune::className<Basis>()), includes);
+
+    auto [lv, isNew] = Dune::Python::insertClass<LocalViewWrapper>(
+        scopedf, "LocalView_" + Dune::className<typename Basis::PrefixPath>(),
+        Dune::Python::GenerateTypeName("Dune::Python::LocalViewWrapper", Dune::MetaType<Basis>()), includes);
+    if (isNew) {
+      lv.def("bind", &LocalViewWrapper::bind);
+      lv.def("unbind", &LocalViewWrapper::unbind);
+      lv.def("index", [](const LocalViewWrapper& localView, int index) { return localView.index(index); });
+      lv.def("__len__", [](LocalViewWrapper& self) -> int { return self.size(); });
+
+      Dune::Python::Functions::registerTree<typename LocalViewWrapper::Tree>(lv);
+      lv.def("tree", [](const LocalViewWrapper& view) { return view.tree(); });
+    }
+  }
+} // namespace Impl
+
+template <class DirichletValues>
+void forwardCorrectFunction(DirichletValues& dirichletValues, const pybind11::function& functor, auto&& cppFunction) {
+  using Basis        = typename DirichletValues::Basis;
+  using Intersection = typename Basis::GridView::Intersection;
+  using BackendType  = typename DirichletValues::BackendType;
+  using MultiIndex   = typename Basis::MultiIndex;
+
+  // Disambiguate by number of arguments
+  pybind11::module inspect_module = pybind11::module::import("inspect");
+  pybind11::object result         = inspect_module.attr("signature")(functor).attr("parameters");
+  size_t numParams                = pybind11::len(result);
+
+  if (numParams == 2) {
+    auto function = functor.template cast<const Impl::FixBoundaryDOFsWithGlobalIndexFunction>();
+    auto lambda   = [&](BackendType& vec, const MultiIndex& indexGlobal) { function(vec.vector(), indexGlobal); };
+    cppFunction(lambda);
+
+  } else if (numParams == 3) {
+    auto lambda = [&](BackendType& vec, int localIndex, auto&& lv) {
+      using SubSpaceBasis = typename std::remove_cvref_t<decltype(lv)>::GlobalBasis;
+      Impl::registerSubSpaceLocalView<SubSpaceBasis>();
+
+      using SubSpaceLocalViewWrapper = Dune::Python::LocalViewWrapper<SubSpaceBasis>;
+      auto lvWrapper                 = SubSpaceLocalViewWrapper(lv);
+
+      auto function =
+          functor.template cast<const Impl::FixBoundaryDOFsWithLocalViewFunction<SubSpaceLocalViewWrapper>>();
+      function(vec.vector(), localIndex, lvWrapper);
+    };
+    cppFunction(lambda);
+
+  } else if (numParams == 4) {
+    auto lambda = [&](BackendType& vec, int localIndex, auto&& lv, const Intersection& intersection) {
+      using SubSpaceBasis = typename std::remove_cvref_t<decltype(lv)>::GlobalBasis;
+      Impl::registerSubSpaceLocalView<SubSpaceBasis>();
+
+      using SubSpaceLocalViewWrapper = Dune::Python::LocalViewWrapper<SubSpaceBasis>;
+      auto lvWrapper                 = SubSpaceLocalViewWrapper(lv);
+
+      auto function = functor.template cast<
+          const Impl::FixBoundaryDOFsWithIntersectionFunction<SubSpaceLocalViewWrapper, Intersection>>();
+      function(vec.vector(), localIndex, lvWrapper, intersection);
+    };
+    cppFunction(lambda);
+
+  } else {
+    DUNE_THROW(Dune::NotImplemented, "fixBoundaryDOFs: A function with this signature is not supported");
+  }
+}
+
 /**
  * \brief Register Python bindings for a DirichletValues class.
  *
  * This function registers Python bindings for a DirichletValues class, allowing it to be used in Python scripts.
  * The registered class will have an initializer that takes a `Basis` object. It exposes several member functions to
  * Python:
- *   - `fixBoundaryDOFs(f)`: Fixes boundary degrees of freedom using a user-defined function `f`.
- *   - `fixBoundaryDOFsUsingLocalView(f)`: Fixes boundary degrees of freedom using a user-defined function `f` with a
- * `LocalView` argument.
- *   - `fixBoundaryDOFsUsingLocalViewAndIntersection(f)`: Fixes boundary degrees of freedom using a user-defined
- * function `f` with `LocalView` and `Intersection` arguments.
- *   - `fixDOFs(f)`: Fixes boundary degrees of freedom using a user-defined function `f` with the boolean vector and
- * the basis as arguments.
+ *  - `fixBoundaryDOFs(f)`: Fixes boundary degrees of freedom using a user-defined function `f` than can be called
+ * with the following arguments:
+ *   -  with the boolean vector and the global index.
+ *   -  with the boolean vector, the local index and the `LocalView`.
+ *   -  with the boolean vector, the local index, the `LocalView` and the `Intersection`.
+ *  - `fixDOFs(f)`: Fixes boundary degrees of freedom using a user-defined function `f` with the basis and the boolean
+ * vector as arguments.
+ *  - `setSingleDOF(i, flag: bool): Fixes or unfixes DOF with index i
+ *  - `isConstrained(i)`: Checks whether index i is constrained
+ *  - `reset()`: Resets the whole container
+ *
+ * The following properties can be accessed:
+ *  - `container`: the underlying container of dirichlet flags (as a const reference)
+ *  - `size`: the size of the underlying basis
+ *  - `fixedDOFsize`: the amount of DOFs currently fixed
  *
  * \tparam DirichletValues The DirichletValues class to be registered.
  * \tparam options Variadic template parameters for additional options when defining the Python class.
@@ -57,60 +156,38 @@ void registerDirichletValues(pybind11::handle scope, pybind11::class_<DirichletV
   using Intersection = typename Basis::GridView::Intersection;
 
   pybind11::module scopedf = pybind11::module::import("dune.functions");
-  typedef Dune::Python::LocalViewWrapper<Basis> LocalViewWrapper;
-  auto includes = Dune::Python::IncludeFiles{"dune/python/functions/globalbasis.hh"};
-  auto lv       = Dune::Python::insertClass<LocalViewWrapper>(
-                scopedf, "LocalView",
-                Dune::Python::GenerateTypeName("Dune::Python::LocalViewWrapper", Dune::MetaType<Basis>()), includes)
-                .first;
+  using LocalViewWrapper   = Dune::Python::LocalViewWrapper<Basis>;
 
-  cls.def(pybind11::init([](const Basis& basis) { return new DirichletValues(basis); }), pybind11::keep_alive<1, 2>());
+  auto includes    = Dune::Python::IncludeFiles{"dune/python/functions/globalbasis.hh"};
+  auto [lv, isNew] = Dune::Python::insertClass<LocalViewWrapper>(
+      scopedf, "LocalView", Dune::Python::GenerateTypeName("Dune::Python::LocalViewWrapper", Dune::MetaType<Basis>()),
+      includes);
 
-  // Eigen::Ref needed due to https://pybind11.readthedocs.io/en/stable/advanced/cast/eigen.html#pass-by-reference
-  cls.def("fixBoundaryDOFs",
-          [](DirichletValues& self, const std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int)>& f) {
-            auto lambda = [&](BackendType& vec, const MultiIndex& indexGlobal) {
-              // we explicitly only allow flat indices
-              f(vec.vector(), indexGlobal[0]);
-            };
-            self.fixBoundaryDOFs(lambda);
-          });
+  if (isNew) {
+    lv.def("bind", &LocalViewWrapper::bind);
+    lv.def("unbind", &LocalViewWrapper::unbind);
+    lv.def("index", [](const LocalViewWrapper& localView, int index) { return localView.index(index); });
+    lv.def("__len__", [](LocalViewWrapper& self) -> int { return self.size(); });
 
-  cls.def("fixBoundaryDOFsUsingLocalView",
-          [](DirichletValues& self,
-             const std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int, LocalViewWrapper&)>& f) {
-            auto lambda = [&](BackendType& vec, int localIndex, LocalView& lv) {
-              auto lvWrapper = LocalViewWrapper(lv.globalBasis());
-              // this can be simplified when
-              // https://gitlab.dune-project.org/staging/dune-functions/-/merge_requests/418 becomes available
-              pybind11::object obj = pybind11::cast(lv.element());
-              lvWrapper.bind(obj);
-              f(vec.vector(), localIndex, lvWrapper);
-            };
-            self.fixBoundaryDOFs(lambda);
-          });
+    Dune::Python::Functions::registerTree<typename LocalViewWrapper::Tree>(lv);
+    lv.def("tree", [](const LocalViewWrapper& view) { return view.tree(); });
+  }
+
+  cls.def(pybind11::init([](const Basis& basis) { return new DirichletValues(basis); }), pybind11::keep_alive<1, 2>());
 
-  cls.def(
-      "fixBoundaryDOFsUsingLocalViewAndIntersection",
-      [](DirichletValues& self,
-         const std::function<void(Eigen::Ref<Eigen::VectorX<bool>>, int, LocalViewWrapper&, const Intersection&)>& f) {
-        auto lambda = [&](BackendType& vec, int localIndex, LocalView& lv, const Intersection& intersection) {
-          auto lvWrapper = LocalViewWrapper(lv.globalBasis());
-          // this can be simplified when
-          // https://gitlab.dune-project.org/staging/dune-functions/-/merge_requests/418 becomes available
-          pybind11::object obj = pybind11::cast(lv.element());
-          lvWrapper.bind(obj);
-          f(vec.vector(), localIndex, lvWrapper, intersection);
-        };
-        self.fixBoundaryDOFs(lambda);
-      });
+  cls.def_property_readonly("container", &DirichletValues::container);
+  cls.def_property_readonly("size", &DirichletValues::size);
+  cls.def("__len__", &DirichletValues::size);
+  cls.def_property_readonly("fixedDOFsize", &DirichletValues::fixedDOFsize);
+  cls.def("isConstrained", [](DirichletValues& self, std::size_t i) -> bool { return self.isConstrained(i); });
+  cls.def("setSingleDOF", [](DirichletValues& self, std::size_t i, bool flag) { self.setSingleDOF(i, flag); });
+  cls.def("isConstrained", [](DirichletValues& self, MultiIndex i) -> bool { return self.isConstrained(i); });
+  cls.def("setSingleDOF", [](DirichletValues& self, MultiIndex i, bool flag) { self.setSingleDOF(i, flag); });
+  cls.def("reset", &DirichletValues::reset);
 
   cls.def("fixDOFs",
           [](DirichletValues& self, const std::function<void(const Basis&, Eigen::Ref<Eigen::VectorX<bool>>)>& f) {
-            auto lambda = [&](const Basis& basis, BackendType& vec) {
-              // we explicitly only allow flat indices
-              f(basis, vec.vector());
-            };
+            auto lambda = [&](const Basis& basis, BackendType& vec) { f(basis, vec.vector()); };
             self.fixDOFs(lambda);
           });
 }
diff --git a/ikarus/python/finiteelements/fe.hh b/ikarus/python/finiteelements/fe.hh
index ec29ad809..34eb02e49 100644
--- a/ikarus/python/finiteelements/fe.hh
+++ b/ikarus/python/finiteelements/fe.hh
@@ -48,7 +48,7 @@ namespace Ikarus::Python {
  * \throws Dune::NotImplemented If the specified resultType is not supported by the finite element.
  */
 template <class FE, class... options>
-void registerCalculateAt(pybind11::handle scope, pybind11::class_<FE, options...> cls, auto restultTypesTuple) {
+void registerCalculateAt(pybind11::handle scope, pybind11::class_<FE, options...> cls, auto resultTypesTuple) {
   using Traits         = typename FE::Traits;
   using FERequirements = typename FE::Requirement;
   cls.def(
@@ -57,7 +57,7 @@ void registerCalculateAt(pybind11::handle scope, pybind11::class_<FE, options...
           std::string resType) {
         Eigen::VectorXd result;
         bool success = false;
-        Dune::Hybrid::forEach(restultTypesTuple, [&]<typename RT>(RT i) {
+        Dune::Hybrid::forEach(resultTypesTuple, [&]<typename RT>(RT i) {
           if (resType == toString(i)) {
             success = true;
             result  = self.template calculateAt<RT::template Rebind>(req, local).asVec();
@@ -113,21 +113,22 @@ void registerFE(pybind11::handle scope, pybind11::class_<FE, options...> cls) {
       pybind11::arg("Requirement"), pybind11::arg("MatrixAffordance"), pybind11::arg("elementMatrix").noconvert());
 
   pybind11::module scopedf = pybind11::module::import("dune.functions");
+  using LocalViewWrapper   = Dune::Python::LocalViewWrapper<FlatBasis>;
 
-  typedef Dune::Python::LocalViewWrapper<FlatBasis> LocalViewWrapper;
-  auto includes = Dune::Python::IncludeFiles{"dune/python/functions/globalbasis.hh"};
-  auto lv       = Dune::Python::insertClass<LocalViewWrapper>(
-                scopedf, "LocalViewWrapper",
-                Dune::Python::GenerateTypeName("Dune::Python::LocalViewWrapperWrapper", Dune::MetaType<FlatBasis>()),
-                includes)
-                .first;
-  lv.def("bind", &LocalViewWrapper::bind);
-  lv.def("unbind", &LocalViewWrapper::unbind);
-  lv.def("index", [](const LocalViewWrapper& localView, int index) { return localView.index(index); });
-  lv.def("__len__", [](LocalViewWrapper& self) -> int { return self.size(); });
-
-  Dune::Python::Functions::registerTree<typename LocalViewWrapper::Tree>(lv);
-  lv.def("tree", [](const LocalViewWrapper& view) { return view.tree(); });
+  auto includes    = Dune::Python::IncludeFiles{"dune/python/functions/globalbasis.hh"};
+  auto [lv, isNew] = Dune::Python::insertClass<LocalViewWrapper>(
+      scopedf, "LocalViewWrapper",
+      Dune::Python::GenerateTypeName("Dune::Python::LocalViewWrapperWrapper", Dune::MetaType<FlatBasis>()), includes);
+
+  if (isNew) {
+    lv.def("bind", &LocalViewWrapper::bind);
+    lv.def("unbind", &LocalViewWrapper::unbind);
+    lv.def("index", [](const LocalViewWrapper& localView, int index) { return localView.index(index); });
+    lv.def("__len__", [](LocalViewWrapper& self) -> int { return self.size(); });
+
+    Dune::Python::Functions::registerTree<typename LocalViewWrapper::Tree>(lv);
+    lv.def("tree", [](const LocalViewWrapper& view) { return view.tree(); });
+  }
 
   cls.def(
       "localView",
diff --git a/ikarus/python/finiteelements/registerferequirements.hh b/ikarus/python/finiteelements/registerferequirements.hh
index 77db64de7..83f2158ca 100644
--- a/ikarus/python/finiteelements/registerferequirements.hh
+++ b/ikarus/python/finiteelements/registerferequirements.hh
@@ -3,6 +3,7 @@
 
 #pragma once
 
+#include <dune/python/common/typeregistry.hh>
 #include <dune/python/pybind11/pybind11.h>
 
 #include <ikarus/python/finiteelements/scalarwrapper.hh>
@@ -22,26 +23,27 @@ void registerFERequirement(pybind11::handle scope, pybind11::class_<FE, options.
       "createRequirement", [](pybind11::object /* self */) { return FERequirements(); },
       pybind11::return_value_policy::copy);
 
-  auto includes              = Dune::Python::IncludeFiles{"ikarus/finiteelements/ferequirements.hh"};
-  auto [lv, isNotRegistered] = Dune::Python::insertClass<FERequirements>(
+  auto includes     = Dune::Python::IncludeFiles{"ikarus/finiteelements/ferequirements.hh"};
+  auto [req, isNew] = Dune::Python::insertClass<FERequirements>(
       scope, "FERequirements", Dune::Python::GenerateTypeName(Dune::className<FERequirements>()), includes);
-  if (isNotRegistered) {
-    lv.def(pybind11::init());
-    lv.def(pybind11::init<SolutionVectorType&, ParameterType&>());
 
-    lv.def(
+  if (isNew) {
+    req.def(pybind11::init());
+    req.def(pybind11::init<SolutionVectorType&, ParameterType&>());
+
+    req.def(
         "insertGlobalSolution",
         [](FERequirements& self, SolutionVectorType solVec) { self.insertGlobalSolution(solVec); },
         "solutionVector"_a.noconvert());
-    lv.def(
+    req.def(
         "globalSolution", [](FERequirements& self) { return self.globalSolution(); },
         pybind11::return_value_policy::reference_internal);
-    lv.def(
+    req.def(
         "insertParameter",
         [](FERequirements& self, ScalarWrapper<double>& parVal) { self.insertParameter(parVal.value()); },
         pybind11::keep_alive<1, 2>(), "parameterValue"_a.noconvert());
 
-    lv.def("parameter", [](const FERequirements& self) { return self.parameter(); });
+    req.def("parameter", [](const FERequirements& self) { return self.parameter(); });
   }
 }
 } // namespace Ikarus::Python
\ No newline at end of file
diff --git a/ikarus/python/test/CMakeLists.txt b/ikarus/python/test/CMakeLists.txt
index 1c869f9d1..ddcf6e290 100644
--- a/ikarus/python/test/CMakeLists.txt
+++ b/ikarus/python/test/CMakeLists.txt
@@ -34,6 +34,17 @@ dune_python_add_test(
   python
 )
 
+dune_python_add_test(
+  NAME
+  pydirichletvalues
+  SCRIPT
+  dirichletvaluetest.py
+  WORKING_DIRECTORY
+  ${CMAKE_CURRENT_SOURCE_DIR}
+  LABELS
+  python
+)
+
 if(HAVE_DUNE_IGA)
   dune_python_add_test(
     NAME
diff --git a/ikarus/python/test/dirichletvaluetest.py b/ikarus/python/test/dirichletvaluetest.py
new file mode 100644
index 000000000..3ab49a59a
--- /dev/null
+++ b/ikarus/python/test/dirichletvaluetest.py
@@ -0,0 +1,112 @@
+# SPDX-FileCopyrightText: 2021-2024 The Ikarus Developers mueller@ibb.uni-stuttgart.de
+# SPDX-License-Identifier: LGPL-3.0-or-later
+
+import debug_info
+
+debug_info.setDebugFlags()
+
+import ikarus as iks
+import dune.grid
+import dune.functions
+
+import math
+import numpy as np
+
+
+def makeGrid():
+    lowerLeft = []
+    upperRight = []
+    elements = []
+    for _ in range(2):
+        lowerLeft.append(-1)
+        upperRight.append(1)
+        elements.append(3)
+
+    return dune.grid.structuredGrid(lowerLeft, upperRight, elements)
+
+
+def testDirichletValues():
+    grid = makeGrid()
+    basis = iks.basis(grid, dune.functions.Power(dune.functions.Lagrange(order=1), 2))
+    basis = basis.flat()
+
+    dirichletValues = iks.dirichletValues(basis)
+
+    assert basis.size() == dirichletValues.size
+    assert basis.size() == len(dirichletValues)
+
+    def fixOneIndex(vec, globalIndex):
+        if globalIndex == 0:
+            vec[globalIndex] = True
+
+    dirichletValues.fixBoundaryDOFs(fixOneIndex)
+
+    # This is equivalent to values[0] == True, but this syntax is discouraged by PEP
+    assert dirichletValues.container[0]
+
+    # Note that the result of localView.index(localIndex) is a multiIndex even for a flat basis, the localIndex is an int
+    def fixAnotherIndexWithLocalView(vec, localIndex, localView):
+        if localView.index(localIndex) == [4]:
+            vec[localView.index(localIndex)] = True
+
+        assert isinstance(localIndex, int)
+        assert isinstance(localView.index(localIndex), list)
+        assert isinstance(localView.index(localIndex)[0], int)
+
+    dirichletValues.fixBoundaryDOFs(fixAnotherIndexWithLocalView)
+    container = dirichletValues.container
+
+    assert sum(container) == 2
+    assert dirichletValues.fixedDOFsize == 2
+
+    def fixTopSide(vec, localIndex, localView, intersection):
+        if intersection.geometry.center[0] == 1.0:
+            vec[localView.index(localIndex)] = True
+
+        assert isinstance(localIndex, int)
+        assert isinstance(localView.index(localIndex), list)
+        assert isinstance(localView.index(localIndex)[0], int)
+
+    dirichletValues.fixBoundaryDOFs(fixTopSide)
+
+    # This assmues a structured grid
+    indicesPerDirection: int = (math.sqrt(grid.size(0)) + 1) * 2
+    assert dirichletValues.fixedDOFsize == 2 + indicesPerDirection
+
+    # This checks whether container is a reference
+    assert sum(container) == dirichletValues.fixedDOFsize
+
+    # Test Subbasis
+    dirichletValues2 = iks.dirichletValues(basis)
+
+    dirichletValues2.fixBoundaryDOFs(fixOneIndex, 0)
+    dirichletValues2.fixBoundaryDOFs(fixAnotherIndexWithLocalView, 0)
+    dirichletValues2.fixBoundaryDOFs(fixTopSide, 0)
+
+    assert dirichletValues2.fixedDOFsize == int(2 + indicesPerDirection / 2)
+
+    dirichletValues2.fixBoundaryDOFs(fixTopSide, 1)
+    assert dirichletValues2.fixedDOFsize == 2 + indicesPerDirection
+
+    dirichletValues2.setSingleDOF(1, True)
+    assert dirichletValues2.fixedDOFsize == 2 + indicesPerDirection + 1
+    assert dirichletValues2.container[1]
+    assert dirichletValues2.isConstrained(1)
+
+    dirichletValues2.setSingleDOF((1), False)  # via MultiIndex
+    assert dirichletValues2.fixedDOFsize == 2 + indicesPerDirection
+    assert not dirichletValues2.isConstrained((1))  # via MultiIndex
+
+    dirichletValues2.reset()
+    assert dirichletValues2.fixedDOFsize == 0
+    assert sum(dirichletValues2.container) == 0
+
+    def fixDOFFunction(basis, vec):
+        vec[:] = np.ones(dirichletValues2.size)
+
+    dirichletValues2.fixDOFs(fixDOFFunction)
+    assert dirichletValues2.fixedDOFsize == dirichletValues2.size
+
+
+if __name__ == "__main__":
+    testDirichletValues()
diff --git a/ikarus/python/test/kltest.py b/ikarus/python/test/kltest.py
index 64abc2c73..e6957c480 100644
--- a/ikarus/python/test/kltest.py
+++ b/ikarus/python/test/kltest.py
@@ -82,7 +82,7 @@ def fixLeftAndRightEdge(vec, localIndex, localView, intersection):
         ):
             vec[localView.index(localIndex)] = True
 
-    dirichletValues.fixBoundaryDOFsUsingLocalViewAndIntersection(fixLeftAndRightEdge)
+    dirichletValues.fixBoundaryDOFs(fixLeftAndRightEdge)
 
     assembler = iks.assembler.sparseFlatAssembler(fes, dirichletValues)
 
diff --git a/ikarus/python/test/linearelastictest.py b/ikarus/python/test/linearelastictest.py
index 5afb65dd0..70e13e9c7 100644
--- a/ikarus/python/test/linearelastictest.py
+++ b/ikarus/python/test/linearelastictest.py
@@ -1,8 +1,8 @@
 # SPDX-FileCopyrightText: 2021-2024 The Ikarus Developers mueller@ibb.uni-stuttgart.de
 # SPDX-License-Identifier: LGPL-3.0-or-later
 
-
 import debug_info
+
 debug_info.setDebugFlags()
 
 import ikarus as iks
@@ -57,20 +57,26 @@ def loadTopEdgePredicate(x):
 
     indexSet = grid.indexSet
     for v in grid.vertices:
-        neumannVertices[indexSet.index(v)]=loadTopEdgePredicate(v.geometry.center)
+        neumannVertices[indexSet.index(v)] = loadTopEdgePredicate(v.geometry.center)
 
     boundaryPatch = iks.utils.boundaryPatch(grid, neumannVertices)
-    nBLoad= iks.finite_elements.neumannBoundaryLoad(boundaryPatch,neumannLoad)
+
+    nBLoad = iks.finite_elements.neumannBoundaryLoad(boundaryPatch, neumannLoad)
 
     linElastic = iks.finite_elements.linearElastic(youngs_modulus=1000, nu=0.2)
-    easF= iks.finite_elements.eas(4)
+    easF = iks.finite_elements.eas(4)
 
     for e in grid.elements:
         if easBool:
-
-            fes.append(iks.finite_elements.makeFE(basisLagrange1,linElastic,easF,vLoad,nBLoad))
+            fes.append(
+                iks.finite_elements.makeFE(
+                    basisLagrange1, linElastic, easF, vLoad, nBLoad
+                )
+            )
         else:
-            fes.append(iks.finite_elements.makeFE(basisLagrange1,linElastic,vLoad,nBLoad))
+            fes.append(
+                iks.finite_elements.makeFE(basisLagrange1, linElastic, vLoad, nBLoad)
+            )
         fes[-1].bind(e)
 
     req = fes[0].createRequirement()
@@ -90,8 +96,8 @@ def loadTopEdgePredicate(x):
 
     forces = np.zeros(8)
     stiffness = np.zeros((8, 8))
-    fes[0].calculateVector(req,iks.VectorAffordance.forces, forces)
-    fes[0].calculateMatrix(req,iks.MatrixAffordance.stiffness, stiffness)
+    fes[0].calculateVector(req, iks.VectorAffordance.forces, forces)
+    fes[0].calculateMatrix(req, iks.MatrixAffordance.stiffness, stiffness)
     fes[0].localView()
 
     dirichletValues = iks.dirichletValues(flatBasis)
@@ -108,8 +114,8 @@ def fixLeftHandEdge(vec, localIndex, localView, intersection):
             vec[localView.index(localIndex)] = True
 
     dirichletValues.fixBoundaryDOFs(fixFirstIndex)
-    dirichletValues.fixBoundaryDOFsUsingLocalView(fixAnotherVertex)
-    dirichletValues.fixBoundaryDOFsUsingLocalViewAndIntersection(fixLeftHandEdge)
+    dirichletValues.fixBoundaryDOFs(fixAnotherVertex)
+    dirichletValues.fixBoundaryDOFs(fixLeftHandEdge)
 
     assembler = iks.assembler.sparseFlatAssembler(fes, dirichletValues)
     assemblerDense = iks.assembler.denseFlatAssembler(fes, dirichletValues)
@@ -134,17 +140,15 @@ def fixLeftHandEdge(vec, localIndex, localView, intersection):
     )
     # Writing results into vtk file
 
-    fileName = "resultdisplacement"+ ("EAS" if easBool else "")
+    fileName = "resultdisplacement" + ("EAS" if easBool else "")
 
-    writer = vtkWriter(
-        grid, fileName, pointData={("displacement", (0, 1)): fx}
-    )
+    vtkWriter(grid, fileName, pointData={("displacement", (0, 1)): fx})
 
     writer2 = vtkUnstructuredGridWriter(grid)
     writer2.addCellData(stressFuncScalar, name="stress")
     writer2.addCellData(stressFuncVec, name="stress2")
 
-    writer2.write(name= "result"+ ("EAS" if easBool else ""))
+    writer2.write(name="result" + ("EAS" if easBool else ""))
 
     # Querying for a different ResultType should result in a runtime error
     try:
diff --git a/ikarus/python/test/nonlinearelastictest.py b/ikarus/python/test/nonlinearelastictest.py
index cc30e5075..8d8f44d7c 100644
--- a/ikarus/python/test/nonlinearelastictest.py
+++ b/ikarus/python/test/nonlinearelastictest.py
@@ -2,6 +2,7 @@
 # SPDX-License-Identifier: LGPL-3.0-or-later
 
 import debug_info
+
 debug_info.setDebugFlags()
 
 import ikarus as iks
@@ -47,11 +48,11 @@ def loadTopEdgePredicate(x):
 
     indexSet = grid.indexSet
     for v in grid.vertices:
-        neumannVertices[indexSet.index(v)]=loadTopEdgePredicate(v.geometry.center)
+        neumannVertices[indexSet.index(v)] = loadTopEdgePredicate(v.geometry.center)
 
     boundaryPatch = iks.utils.boundaryPatch(grid, neumannVertices)
 
-    nBLoad= iks.finite_elements.neumannBoundaryLoad(boundaryPatch,neumannLoad)
+    nBLoad = iks.finite_elements.neumannBoundaryLoad(boundaryPatch, neumannLoad)
 
     svk = iks.materials.StVenantKirchhoff(E=1000, nu=0.3)
 
@@ -61,7 +62,9 @@ def loadTopEdgePredicate(x):
 
     fes = []
     for e in grid.elements:
-        fes.append(iks.finite_elements.makeFE(basisLagrange1,nonLinElastic,vLoad,nBLoad))
+        fes.append(
+            iks.finite_elements.makeFE(basisLagrange1, nonLinElastic, vLoad, nBLoad)
+        )
         fes[-1].bind(e)
 
     dirichletValues = iks.dirichletValues(flatBasis)
@@ -78,8 +81,8 @@ def fixLeftHandEdge(vec, localIndex, localView, intersection):
             vec[localView.index(localIndex)] = True
 
     dirichletValues.fixBoundaryDOFs(fixFirstIndex)
-    dirichletValues.fixBoundaryDOFsUsingLocalView(fixAnotherVertex)
-    dirichletValues.fixBoundaryDOFsUsingLocalViewAndIntersection(fixLeftHandEdge)
+    dirichletValues.fixBoundaryDOFs(fixAnotherVertex)
+    dirichletValues.fixBoundaryDOFs(fixLeftHandEdge)
 
     assembler = iks.assembler.sparseFlatAssembler(fes, dirichletValues)
 
@@ -91,13 +94,11 @@ def fixLeftHandEdge(vec, localIndex, localView, intersection):
 
     def energy(dRedInput):
         feReq = fes[0].createRequirement()
-        feReq.insertParameter( lambdaLoad)
+        feReq.insertParameter(lambdaLoad)
         dBig = assembler.createFullVector(dRedInput)
-        feReq.insertGlobalSolution( dBig)
+        feReq.insertGlobalSolution(dBig)
         feReq.globalSolution()
-        return assembler.scalar(
-            feReq, iks.ScalarAffordance.mechanicalPotentialEnergy
-        )
+        return assembler.scalar(feReq, iks.ScalarAffordance.mechanicalPotentialEnergy)
 
     def gradient(dRedInput):
         feReq = fes[0].createRequirement()
@@ -105,7 +106,8 @@ def gradient(dRedInput):
         dBig = assembler.createFullVector(dRedInput)
         feReq.insertGlobalSolution(dBig)
         return assembler.vector(
-            feReq, iks.VectorAffordance.forces, iks.DBCOption.Reduced)
+            feReq, iks.VectorAffordance.forces, iks.DBCOption.Reduced
+        )
 
     def hess(dRedInput):
         feReq = fes[0].createRequirement()
@@ -114,7 +116,7 @@ def hess(dRedInput):
         feReq.insertGlobalSolution(dBig)
         return assembler.matrix(
             feReq, iks.MatrixAffordance.stiffness, iks.DBCOption.Reduced
-        ).todense() # this is slow, but for this test we don't care
+        ).todense()  # this is slow, but for this test we don't care
 
     resultd = minimize(energy, x0=dRed, options={"disp": True}, tol=1e-14)
     resultd2 = minimize(
@@ -137,6 +139,6 @@ def hess(dRedInput):
 
     feReq = fes[0].createRequirement()
     fullD = assembler.createFullVector(resultd2.x)
-    feReq.insertGlobalSolution( fullD)
+    feReq.insertGlobalSolution(fullD)
 
     res1 = fes[0].calculateAt(feReq, np.array([0.5, 0.5]), "PK2Stress")
diff --git a/ikarus/python/test/testmaterials.py b/ikarus/python/test/testmaterials.py
index 82ecfdae1..1b150ec30 100644
--- a/ikarus/python/test/testmaterials.py
+++ b/ikarus/python/test/testmaterials.py
@@ -2,6 +2,7 @@
 # SPDX-License-Identifier: LGPL-3.0-or-later
 
 import debug_info
+
 debug_info.setDebugFlags()
 
 import ikarus as iks
diff --git a/ikarus/utils/dirichletvalues.hh b/ikarus/utils/dirichletvalues.hh
index 18ac9da4e..c4608e407 100644
--- a/ikarus/utils/dirichletvalues.hh
+++ b/ikarus/utils/dirichletvalues.hh
@@ -18,6 +18,8 @@
 
 #include <dune/functions/backends/istlvectorbackend.hh>
 #include <dune/functions/functionspacebases/boundarydofs.hh>
+#include <dune/functions/functionspacebases/flatmultiindex.hh>
+#include <dune/functions/functionspacebases/subspacebasis.hh>
 
 #include <Eigen/Core>
 
@@ -84,23 +86,29 @@ public:
    * degrees of freedom and the usual arguments of `Dune::Functions::forEachBoundaryDOF`.
    *
    * \param f A callback function
+   * \param treePath An optional argument specifying a tree path to a subspacebasis, e.g. Dune::Indices::_0
    */
-  template <typename F>
-  void fixBoundaryDOFs(F&& f) {
+  template <typename F, typename TreePath = Dune::TypeTree::HybridTreePath<>>
+  void fixBoundaryDOFs(F&& f, TreePath&& treePath = {}) {
+    using namespace Dune::Functions;
+    using SubSpaceLocalView =
+        typename std::remove_cvref_t<decltype(subspaceBasis(basis_, std::forward<TreePath>(treePath)))>::LocalView;
+
     if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, typename Basis::MultiIndex>) {
       auto lambda = [&](auto&& indexGlobal) { f(dirichletFlagsBackend_, indexGlobal); };
-      Dune::Functions::forEachBoundaryDOF(basis_, lambda);
-    } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, typename Basis::LocalView>) {
+      Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
+    } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, SubSpaceLocalView>) {
       auto lambda = [&](auto&& localIndex, auto&& localView) { f(dirichletFlagsBackend_, localIndex, localView); };
-      Dune::Functions::forEachBoundaryDOF(basis_, lambda);
-    } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, typename Basis::LocalView,
+      Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
+    } else if constexpr (Concepts::IsFunctorWithArgs<F, BackendType, int, SubSpaceLocalView,
                                                      typename Basis::GridView::Intersection>) {
       auto lambda = [&](auto&& localIndex, auto&& localView, auto&& intersection) {
         f(dirichletFlagsBackend_, localIndex, localView, intersection);
       };
-      Dune::Functions::forEachBoundaryDOF(basis_, lambda);
-    } else
-      DUNE_THROW(Dune::IOError, "Invalid callback function passed to fixBoundaryDOFs");
+      Dune::Functions::forEachBoundaryDOF(subspaceBasis(basis_, std::forward<TreePath>(treePath)), lambda);
+    } else {
+      static_assert(Dune::AlwaysFalse<F>(), "fixBoundaryDOFs: A function with this signature is not supported");
+    }
   }
 
   /**
@@ -117,11 +125,37 @@ public:
   }
 
   /**
-   * \brief Function to fix (set boolean values to true or false) of degrees of freedom
+   * \brief Fixes and unfixes (set boolean value to true or false) a specific degree of freedom
+   *
+   * \param i An index indicating the DOF number to be fixed
+   * \param flag Boolean indicating whether the DOF should fixed or not
+   */
+  template <typename MultiIndex>
+  requires(not std::integral<MultiIndex>)
+  void setSingleDOF(const MultiIndex i, bool flag) {
+    dirichletFlagsBackend_[i] = flag;
+  }
+
+  /**
+   * \brief Fixes or unfixes (set boolean value to true or false) a specific degree of freedom
    *
    * \param i An index indicating the DOF number to be fixed
+   * \param flag Boolean indicating whether the DOF should fixed or not
    */
-  void fixIthDOF(typename Basis::MultiIndex i) { dirichletFlagsBackend_[i] = true; }
+
+  void setSingleDOF(std::size_t i, bool flag)
+  requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
+  {
+    dirichletFlags_[i] = flag;
+  }
+
+  /**
+   * \brief Resets all degrees of freedom
+   */
+  void reset() {
+    std::fill(dirichletFlags_.begin(), dirichletFlags_.end(), false);
+    inhomogeneousBoundaryVectorDummy_.setZero(static_cast<Eigen::Index>(basis_.size()));
+  }
 
   /* \brief Returns the local basis object */
   const auto& basis() const { return basis_; }
@@ -134,7 +168,11 @@ public:
   }
 
   /* \brief Returns a boolean values, if the i-th degree of freedom is constrained */
-  [[nodiscard]] bool isConstrained(std::size_t i) const { return dirichletFlags_[i]; }
+  [[nodiscard]] bool isConstrained(std::size_t i) const
+  requires(std::same_as<typename Basis::MultiIndex, Dune::Functions::FlatMultiIndex<size_t>>)
+  {
+    return dirichletFlags_[i];
+  }
 
   /* \brief Returns how many degrees of freedoms are fixed */
   auto fixedDOFsize() const { return std::ranges::count(dirichletFlags_, true); }
diff --git a/python/ikarus/dirichlet_values.py b/python/ikarus/dirichlet_values.py
index ee61fa5f1..3a2e8c063 100644
--- a/python/ikarus/dirichlet_values.py
+++ b/python/ikarus/dirichlet_values.py
@@ -1,9 +1,38 @@
 # SPDX-FileCopyrightText: 2021-2024 The Ikarus Developers mueller@ibb.uni-stuttgart.de
 # SPDX-License-Identifier: LGPL-3.0-or-later
 
+import types
 from dune.common.hashit import hashIt
 from .generator import MySimpleGenerator
 
+from io import StringIO
+from dune.generator.algorithm import run
+
+
+def __fixBoundaryDOFs(dirichletValues):
+    def __fixBoundaryDOFsFunc(dirichletValues, f, *args: int):
+        prefixPathTypeName = "Dune::TypeTree::HybridTreePath<"
+        prefixPathTypeName += ",".join(
+            "Dune::index_constant<" + str(i) + ">" for i in args
+        )
+        prefixPathTypeName += ">"
+
+        runCode = """
+        #include <ikarus/python/dirichletvalues/dirichletvalues.hh> 
+        template <typename DV>
+        void fixBoundaryDofs(DV& dirichletValues, const pybind11::function& functor)
+        {{
+        Ikarus::Python::forwardCorrectFunction(dirichletValues, functor, 
+        [&](auto&& functor_) {{ 
+        return dirichletValues.fixBoundaryDOFs(functor_, {prefixPathType}{{}});
+          }});
+        }}
+        """.format(prefixPathType=prefixPathTypeName)
+        return run("fixBoundaryDofs", StringIO(runCode), dirichletValues, f)
+
+    return __fixBoundaryDOFsFunc
+
+
 def dirichletValues(basis):
     """
     @brief Creates a Dirichlet values handler for the given basis.
@@ -21,6 +50,10 @@ def dirichletValues(basis):
     includes += ["ikarus/python/dirichletvalues/dirichletvalues.hh"]
     moduleName = "dirichletValues_" + hashIt(element_type)
     module = generator.load(
-        includes=includes, typeName=element_type, moduleName=moduleName
+        includes=includes, typeName=element_type, moduleName=moduleName, dynamicAttr=True
     )
-    return module.DirichletValues(basis)
+
+    dv = module.DirichletValues(basis)
+    dv.fixBoundaryDOFs = types.MethodType(__fixBoundaryDOFs(dv), dv)
+
+    return dv
diff --git a/tests/src/testassembler.cpp b/tests/src/testassembler.cpp
index 798a96d7f..c46c8dfb9 100644
--- a/tests/src/testassembler.cpp
+++ b/tests/src/testassembler.cpp
@@ -89,7 +89,7 @@ auto SimpleAssemblersTest(const PreBasis& preBasis) {
     } else {
       const auto fixIndices = utils::globalIndexFromGlobalPosition(basis.flat(), pos);
       for (const auto idx : fixIndices)
-        dirichletValues.fixIthDOF(idx);
+        dirichletValues.setSingleDOF(idx, true);
     }
 
     using SparseAssembler = SparseFlatAssembler<decltype(fes), decltype(dirichletValues)>;
diff --git a/tests/src/testdirichletvalue.cpp b/tests/src/testdirichletvalue.cpp
index 2559dc6ed..1a1ad5984 100644
--- a/tests/src/testdirichletvalue.cpp
+++ b/tests/src/testdirichletvalue.cpp
@@ -2,8 +2,6 @@
 // SPDX-License-Identifier: LGPL-3.0-or-later
 #include <config.h>
 
-#include "testhelpers.hh"
-
 #include <vector>
 
 #include <dune/common/float_cmp.hh>
@@ -29,7 +27,7 @@
 using Dune::TestSuite;
 
 static auto dirichletBCTest() {
-  TestSuite t("dirichletBCTest");
+  TestSuite t("DirichletValueTest");
   using Grid = Dune::YaspGrid<2>;
 
   const double Lx                         = 4.0;
@@ -43,22 +41,120 @@ static auto dirichletBCTest() {
   using namespace Dune::Functions::BasisFactory;
   auto basis = Ikarus::makeBasis(gridView, power<2>(lagrange<1>()));
 
-  auto basisP = std::make_shared<const decltype(basis)>(basis);
+  auto flatBasis = basis.flat();
 
-  Ikarus::DirichletValues dirichletValues1(basisP->flat());
+  Ikarus::DirichletValues dirichletValues1(flatBasis);
   dirichletValues1.fixDOFs([](auto& basis_, auto& dirichFlags) {
     Dune::Functions::forEachBoundaryDOF(basis_, [&](auto&& indexGlobal) { dirichFlags[indexGlobal] = true; });
   });
 
-  Ikarus::DirichletValues dirichletValues2(basisP->flat());
+  Ikarus::DirichletValues dirichletValues2(flatBasis);
   dirichletValues2.fixBoundaryDOFs([](auto& dirichFlags, auto&& indexGlobal) { dirichFlags[indexGlobal] = true; });
 
-  for (std::size_t i = 0; i < basisP->flat().size(); ++i)
+  for (std::size_t i = 0; i < flatBasis.size(); ++i)
     t.check(dirichletValues1.isConstrained(i) == dirichletValues2.isConstrained(i))
         << "Different dirichlet value creations (dirichletValues1 and dirichletValues2) didn't provide the same "
            "result. Index: i="
         << i;
 
+  auto testSubSpaceBasis = [&](auto&& tree1, auto&& tree2) {
+    Ikarus::DirichletValues dirichletValues_SSB(flatBasis);
+    dirichletValues_SSB.fixBoundaryDOFs([](auto& dirichFlags, auto&& indexGlobal) { dirichFlags[indexGlobal] = true; },
+                                        tree1);
+
+    t.check(dirichletValues2.fixedDOFsize() == dirichletValues_SSB.fixedDOFsize() * 2)
+        << "DirichletValues with subspace basis should have half as many fixed DOFs as the full basis, but has "
+        << dirichletValues_SSB.fixedDOFsize() << ", where as the full basis has " << dirichletValues2.fixedDOFsize()
+        << " fixed DOFs";
+
+    dirichletValues_SSB.fixBoundaryDOFs([](auto& dirichFlags, auto&& indexGlobal) { dirichFlags[indexGlobal] = true; },
+                                        tree2);
+
+    for (std::size_t i = 0; i < flatBasis.size(); ++i)
+      t.check(dirichletValues_SSB.isConstrained(i) == dirichletValues2.isConstrained(i))
+          << "Different dirichlet value creations with subspace basis didn't provide the same result as with full "
+             "basis. "
+             "Index: i="
+          << i;
+  };
+
+  // Static tree path
+  auto treePath0 = Dune::Indices::_0;
+  auto treePath1 = Dune::Indices::_1;
+  testSubSpaceBasis(treePath0, treePath1);
+
+  // Dynamic tree path
+  auto treePath0d = 0;
+  auto subBasis1d = 1;
+  testSubSpaceBasis(treePath0d, subBasis1d);
+
+  // Test with Intersection
+  auto fixLambda = [](auto& dirichletFlags, auto&& localIndex, auto&& localView, auto&& intersection) {
+    if (intersection.geometry().center()[0] > 4 - 1e-8)
+      dirichletFlags[localView.index(localIndex)] = true;
+  };
+  Ikarus::DirichletValues dirichletValues5(flatBasis);
+  dirichletValues5.fixBoundaryDOFs(fixLambda);
+
+  Ikarus::DirichletValues dirichletValues_SSB2(flatBasis);
+  dirichletValues_SSB2.fixBoundaryDOFs(fixLambda, treePath0);
+
+  t.check(dirichletValues5.fixedDOFsize() == dirichletValues_SSB2.fixedDOFsize() * 2)
+      << "DirichletValues with subspace basis should have half as many fixed DOFs as the full basis, but has "
+      << dirichletValues_SSB2.fixedDOFsize() << ", where as the full basis has " << dirichletValues5.fixedDOFsize()
+      << " fixed DOFs";
+
+  dirichletValues_SSB2.fixBoundaryDOFs(fixLambda, treePath1);
+
+  t.check(dirichletValues_SSB2.fixedDOFsize() > 0);
+  t.check(dirichletValues5.fixedDOFsize() == dirichletValues_SSB2.fixedDOFsize())
+      << "DirichletValues with subspace basis should have as many fixed DOFs as the full basis, but has "
+      << dirichletValues_SSB2.fixedDOFsize() << ", where as the full basis has " << dirichletValues5.fixedDOFsize()
+      << " fixed DOFs";
+
+  for (std::size_t i = 0; i < flatBasis.size(); ++i)
+    t.check(dirichletValues_SSB2.isConstrained(i) == dirichletValues5.isConstrained(i))
+        << "Different dirichlet value creations with subspace basis didn't provide the same result as with full basis. "
+           "Index: i="
+        << i;
+
+  auto sum        = [](const auto& container_) { return std::accumulate(container_.begin(), container_.end(), 0); };
+  auto manual_sum = [](const auto& dv) {
+    int sum = 0;
+    for (auto i : Dune::range(dv.size()))
+      if (dv.isConstrained(i))
+        ++sum;
+    return sum;
+  };
+
+  // Test container
+  auto& container = dirichletValues2.container();
+  static_assert(std::is_reference_v<decltype(container)>, "container should be a reference");
+  static_assert(std::is_const_v<std::remove_reference_t<decltype(container)>>, "container should be const");
+
+  t.check(std::ranges::all_of(container, [](auto&& flag) { return flag; })) << "All values should be fixed";
+
+  auto sumPre    = sum(container);
+  auto sumManual = manual_sum(dirichletValues2);
+
+  t.check(sumManual == sumPre)
+      << "Summing over container should yield the same amount of fixed DOFs then checking manually";
+
+  dirichletValues2.setSingleDOF(decltype(flatBasis)::MultiIndex{1}, false);
+  t.check(sum(container) == sumPre - 1) << "The sum of fixed DOFs should be one less then after unfixing one";
+
+  t.check(manual_sum(dirichletValues2) == sumPre - 1)
+      << "Summing over container should yield the same amount of fixed DOFs then checking manually";
+
+  dirichletValues2.setSingleDOF(2, false);
+  t.check(sum(container) == sumPre - 2) << "The sum of fixed DOFs should be two less then after unfixing two";
+
+  // Reset
+  dirichletValues2.reset();
+  t.check(sum(container) == 0) << "After resetting container should have only false entries";
+  t.check(manual_sum(dirichletValues2) == 0) << "After resetting, all DOFs should be false";
+
+  // Inhomogenious Boundary Conditions
   auto inhomogeneousDisplacement = []<typename T>(const auto& globalCoord, const T& lambda) {
     Eigen::Vector<T, 2> localInhomogeneous;
     if (globalCoord[0] > 4 - 1e-8) {
@@ -97,7 +193,7 @@ static auto dirichletBCTest() {
           << "Values differ dispDerivs[i]: " << dispDerivs[globalIndex0];
     }
   };
-  Dune::Functions::forEachBoundaryDOF(basisP->flat(), lambdaCheck);
+  Dune::Functions::forEachBoundaryDOF(flatBasis, lambdaCheck);
 
   // Check that we can store lambda from python
   std::string inhomogeneousDisplacementFunction =
@@ -132,7 +228,7 @@ static auto dirichletBCTest() {
     return Dune::toEigen(pythonFunc(globalCoord, lambda_));
   };
   //
-  Ikarus::DirichletValues dirichletValues3(basisP->flat());
+  Ikarus::DirichletValues dirichletValues3(flatBasis);
   dirichletValues3.fixBoundaryDOFs([](auto& dirichFlags, auto&& indexGlobal) { dirichFlags[indexGlobal] = true; });
 
   Eigen::VectorXd disps2, dispDerivs2;
@@ -143,7 +239,7 @@ static auto dirichletBCTest() {
   t.check(disps2.isApprox(dispDerivs2 * lambda3));
 
   // Check if all boundary DOFs found manually using obtainLagrangeNodePositions is the same as using forEachBoundaryDOF
-  Ikarus::DirichletValues dirichletValues4(basisP->flat());
+  Ikarus::DirichletValues dirichletValues4(flatBasis);
   constexpr double tol = 1e-8;
   auto localView       = basis.flat().localView();
   for (auto& ele : elements(gridView)) {
@@ -156,11 +252,11 @@ static auto dirichletBCTest() {
           (std::abs(nodalPos[i][1]) < tol) or (std::abs(nodalPos[i][1] - Ly) < tol))
         for (auto fixedDirection = 0; fixedDirection < 2; ++fixedDirection) {
           auto fixIndex = localView.index(localView.tree().child(fixedDirection).localIndex(i));
-          dirichletValues4.fixIthDOF(fixIndex);
+          dirichletValues4.setSingleDOF(fixIndex, true);
         }
   }
 
-  for (std::size_t i = 0; i < basisP->flat().size(); ++i)
+  for (std::size_t i = 0; i < flatBasis.size(); ++i)
     t.check(dirichletValues1.isConstrained(i) == dirichletValues4.isConstrained(i))
         << "Different dirichlet value creations (dirichletValues1 and dirichletValues4) didn't provide the same "
            "result. Index: i="