diff --git a/README.md b/README.md
index ddec592..363d602 100644
--- a/README.md
+++ b/README.md
@@ -49,6 +49,19 @@ The code is available under the [Apache 2 License](http://www.apache.org/license
 
 If you use the code in your project/product, please drop us a note – we are always interested in learning about new applications!
 
+# Authors & acknowledgements
+
+Authors:
+
+ * Martin Kutz, FU Berlin
+ * [Kaspar Fischer](http://github.com/hbf), ETH Zurich
+ * [Bernd Gärtner](http://www.inf.ethz.ch/personal/gaertner/), ETH Zurich
+
+Many thanks go to the following people who have – sometimes substantially – contributed to the code:
+
+
+ * Thomas Otto (University of Hamburg) for [submitting several compiler fixes](https://github.com/hbf/miniball/issues/3) (g++ 4.7 and 4.5 on SUSE Linux 12.2 and 11.3) and for [introducing generic point and point coordinate accessors](https://github.com/hbf/miniball/pull/5) in the code.
+
 # Links
   * For small dimensions like 2D or 3D, [Bernd Gärtner's code](http://www.inf.ethz.ch/personal/gaertner/miniball.html), which is based on Welzl's algorithm, may be faster.
   * The [Computational Geometry Algorithms Library (CGAL)](http://www.cgal.org/) contains both floating-point and _arbitrary-precision arithmetic_ implementations of several bounding sphere algorithms. Among then, there is an algorithm to compute the minsphere of _spheres_ (not just points as input). See the [Chapter on Bounding Volumes](http://www.cgal.org/Manual/latest/doc_html/cgal_manual/Bounding_volumes/Chapter_main.html).
diff --git a/cpp/main/Seb.C b/cpp/main/Seb.C
index e033478..e234172 100644
--- a/cpp/main/Seb.C
+++ b/cpp/main/Seb.C
@@ -9,10 +9,12 @@
 #include <numeric>
 #include <algorithm>
 
+#include "Seb.h" // Note: header included for better syntax highlighting in some IDEs.
+
 namespace SEB_NAMESPACE {
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::allocate_resources()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::allocate_resources()
   {
     center            = new Float[dim];
     center_to_aff     = new Float[dim];
@@ -20,8 +22,8 @@ namespace SEB_NAMESPACE {
     lambdas           = new Float[dim+1];
   }
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::deallocate_resources()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::deallocate_resources()
   {
     delete[] center;
     delete[] center_to_aff;
@@ -32,8 +34,8 @@ namespace SEB_NAMESPACE {
       delete support;
   }
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::init_ball()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::init_ball()
   // Precondition: |S| > 0
   // Sets up the search ball with an arbitrary point of S as center
   // and with with exactly one of the points farthest from center in
@@ -67,15 +69,15 @@ namespace SEB_NAMESPACE {
     // initialize support to the farthest point:
     if (support != NULL)
       delete support;
-    support = new Subspan<Float>(dim,S,farthest);
+    support = new Subspan<Float, Pt, PointAccessor>(dim,S,farthest);
     
     // statistics:
     // initialize entry-counters to zero:
     SEB_STATS(entry_count = std::vector<int>(S.size(),0));
   }
   
-  template<typename Float>
-  bool Smallest_enclosing_ball<Float>::successful_drop()
+  template<typename Float, class Pt, class PointAccessor>
+  bool Smallest_enclosing_ball<Float, Pt, PointAccessor>::successful_drop()
   // Precondition: center lies in aff(support).
   // If center doesn't already lie in conv(support) and is thus
   // not optimal yet, successful_drop() elects a suitable point k to
@@ -104,8 +106,8 @@ namespace SEB_NAMESPACE {
     return false;
   }
   
-  template<typename Float>
-  Float Smallest_enclosing_ball<Float>::find_stop_fraction(int& stopper)
+  template<typename Float, class Pt, class PointAccessor>
+  Float Smallest_enclosing_ball<Float, Pt, PointAccessor>::find_stop_fraction(int& stopper)
   // Given the center of the current enclosing ball and the
   // walking direction center_to_aff, determine how much we can walk
   // into this direction without losing a point from S.  The (positive)
@@ -163,8 +165,8 @@ namespace SEB_NAMESPACE {
   }
   
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::update()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::update()
   // The main function containing the main loop.
   // Iteratively, we compute the point in support that is closest
   // to the current center and then walk towards this target as far
@@ -274,8 +276,8 @@ namespace SEB_NAMESPACE {
     }
   }
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::verify()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::verify()
   {
     using std::inner_product;
     using std::abs;
@@ -351,8 +353,8 @@ namespace SEB_NAMESPACE {
   }
   
   
-  template<typename Float>
-  void Smallest_enclosing_ball<Float>::test_affine_stuff()
+  template<typename Float, class Pt, class PointAccessor>
+  void Smallest_enclosing_ball<Float, Pt, PointAccessor>::test_affine_stuff()
   {
     using std::cout;
     using std::endl;
@@ -364,7 +366,7 @@ namespace SEB_NAMESPACE {
     cout << S.size() << " points in " << dim << " dimensions" << endl;
     
     if (!is_empty()) {
-      support = new Subspan<Float>(dim,S,0);
+      support = new Subspan<Float,PointAccessor,Pt>(dim,S,0);
       cout << "initializing affine subspace with point S[0]" << endl;
       
       direction = new Float[dim];
@@ -426,8 +428,8 @@ namespace SEB_NAMESPACE {
   }
   
   
-  template<typename Float>
-  const Float Smallest_enclosing_ball<Float>::Eps = Float(1e-14);
+  template<typename Float, class Pt, class PointAccessor>
+  const Float Smallest_enclosing_ball<Float, Pt, PointAccessor>::Eps = Float(1e-14);
   
 } // namespace SEB_NAMESPACE
 
diff --git a/cpp/main/Seb.h b/cpp/main/Seb.h
index 619942c..f2b28c9 100644
--- a/cpp/main/Seb.h
+++ b/cpp/main/Seb.h
@@ -8,32 +8,38 @@
 
 #include <vector>
 #include "Seb_configure.h"
+#include "Seb_point.h"
 #include "Subspan.h"
 
 namespace SEB_NAMESPACE {
+
+  // template arguments:
+  // Float must be a floating point data type for which * / - + are defined
+  // Pt[i] must return the i-th coordinate as a Float
+  // PointAccessor[j] must return the j-th point in the data set as Pt and
+  // size_t size() returns the size of the data set
   
-  template<typename Float>
+  template<typename Float, class Pt = Point<Float>, class PointAccessor = std::vector<Pt> >
   class Smallest_enclosing_ball
   // An instance of class Smallest_enclosing_ball<Float> represents
   // the smallest enclosing ball of a set S of points.  Initially, the
   // set S is empty; you can add points by calling insert().
   {
-  private: // internal representation of points:
-    typedef Point<Float> Pt;
-    
   public: // iterator-type to iterate over the center coordinates of
 	  // the miniball (cf. center_begin() below):
     typedef Float *Coordinate_iterator;
     
   public: // construction and destruction:
-    
-    Smallest_enclosing_ball(unsigned int d)
-    // Constructs an instance representing the miniball of the empty
-    // set S={}.  The dimension of the ambient space is fixed to d for
+
+    Smallest_enclosing_ball(unsigned int d, PointAccessor &P)
+    // Constructs an instance representing the miniball of points from
+    // set S.  The dimension of the ambient space is fixed to d for
     // lifetime of the instance.
-    : dim(d), up_to_date(true), support(NULL)
+    : dim(d), S(P), up_to_date(true), support(NULL)
     {
       allocate_resources();
+      SEB_ASSERT(!is_empty());
+      update();
     }
     
     ~Smallest_enclosing_ball()
@@ -43,12 +49,12 @@ namespace SEB_NAMESPACE {
     
   public: // modification:
     
-    template<typename InputIterator>
-    void insert(InputIterator first)
-    // Inserts the point p into the instance's set S.  The point p is
-    // specified by its d coordinates from the range [first,first+d).
+    void invalidate()
+    // Notifies the instance that the underlying point set S (passed to the constructor
+    // of this instance as parameter P) has changed. This will cause the miniball to
+    // be recomputed lazily (i.e., when you call for example radius(), the recalculation
+    // will be triggered).
     {
-      S.push_back(Pt(dim,first));    //  TODO!  open design decision
       up_to_date = false;
     }
     
@@ -141,12 +147,12 @@ namespace SEB_NAMESPACE {
     
   private: // member fields:
     unsigned int dim;                 // dimension of the amient space
-    std::vector<Pt> S;                // set S of inserted points
+    PointAccessor &S;                 // set S of inserted points
     bool up_to_date;                  // whether the miniball has
-    // already been computed
+                                      // already been computed
     Float *center;                    // center of the miniball
     Float radius_, radius_square;     // squared radius of the miniball
-    Subspan<Float> *support;          // the points that lie on the current
+    Subspan<Float, Pt, PointAccessor> *support;          // the points that lie on the current
     // boundary and "support" the ball;
     // the essential structure for update()
     
diff --git a/cpp/main/Subspan.C b/cpp/main/Subspan.C
index d386d4c..4ebb31b 100644
--- a/cpp/main/Subspan.C
+++ b/cpp/main/Subspan.C
@@ -9,6 +9,8 @@
 #include <numeric>
 #include <cmath>
 
+#include "Subspan.h" // Note: header included for better syntax highlighting in some IDEs.
+
 // The point members[r] is called the origin; we use the following macro
 // to increase the readibilty of the code:
 #define SEB_AFFINE_ORIGIN S[members[r]]
@@ -46,8 +48,8 @@ namespace SEB_NAMESPACE {
     }
   }
   
-  template<typename Float>
-  Subspan<Float>::Subspan(unsigned int dim, const std::vector<Pt>& S, int index)
+  template<typename Float, class Pt, class PointAccessor>
+  Subspan<Float, Pt, PointAccessor>::Subspan(unsigned int dim, const PointAccessor& S, int index)
   : S(S), membership(S.size()), dim(dim), members(dim+1)
   {
     // allocate storage for Q, R, u, and w:
@@ -71,8 +73,8 @@ namespace SEB_NAMESPACE {
     SEB_LOG ("ranks",r << std::endl);
   }
   
-  template<typename Float>
-  Subspan<Float>::~Subspan()
+  template<typename Float, class Pt, class PointAccessor>
+  Subspan<Float, Pt, PointAccessor>::~Subspan()
   {
     for (unsigned int i=0; i<dim; ++i) {
       delete[] Q[i];
@@ -84,8 +86,8 @@ namespace SEB_NAMESPACE {
     delete[] w;
   }
   
-  template<typename Float>
-  void Subspan<Float>::add_point(int index) {
+  template<typename Float, class Pt, class PointAccessor>
+  void Subspan<Float, Pt, PointAccessor>::add_point(int index) {
     SEB_ASSERT(!is_member(index));
     
     // compute S[i] - origin into u:
@@ -105,8 +107,8 @@ namespace SEB_NAMESPACE {
     SEB_LOG ("ranks",r << std::endl);
   }
   
-  template<typename Float>
-  void Subspan<Float>::remove_point(const unsigned int local_index) {
+  template<typename Float, class Pt, class PointAccessor>
+  void Subspan<Float, Pt, PointAccessor>::remove_point(const unsigned int local_index) {
     SEB_ASSERT(is_member(global_index(local_index)) && size() > 1);
     
     membership[global_index(local_index)] = false;
@@ -144,10 +146,10 @@ namespace SEB_NAMESPACE {
     }
   }
   
-  template<typename Float>
+  template<typename Float, class Pt, class PointAccessor>
   template<typename RandomAccessIterator1,
   typename RandomAccessIterator2>
-  Float Subspan<Float>::
+  Float Subspan<Float, Pt, PointAccessor>::
   shortest_vector_to_span(RandomAccessIterator1 p,
                           RandomAccessIterator2 w)
   {
@@ -167,8 +169,8 @@ namespace SEB_NAMESPACE {
     return inner_product(w,w+dim,w,Float(0));
   }
   
-  template<typename Float>
-  Float Subspan<Float>::representation_error()
+  template<typename Float, class Pt, class PointAccessor>
+  Float Subspan<Float, Pt, PointAccessor>::representation_error()
   {
     using std::abs;
     
@@ -199,10 +201,10 @@ namespace SEB_NAMESPACE {
     return max;
   }
   
-  template<typename Float>
+  template<typename Float, class Pt, class PointAccessor>
   template<typename RandomAccessIterator1,
   typename RandomAccessIterator2>
-  void Subspan<Float>::
+  void Subspan<Float, Pt, PointAccessor>::
   find_affine_coefficients(RandomAccessIterator1 p,
                            RandomAccessIterator2 lambdas)
   {
@@ -234,8 +236,8 @@ namespace SEB_NAMESPACE {
     *(lambdas+r) = origin_lambda;
   }
   
-  template<typename Float>
-  void Subspan<Float>::append_column()
+  template<typename Float, class Pt, class PointAccessor>
+  void Subspan<Float, Pt, PointAccessor>::append_column()
   // Appends the new column u (which is a member of this instance) to
   // the right of "A = QR", updating Q and R.  It assumes r to still
   // be the old value, i.e., the index of the column used now for
@@ -274,8 +276,8 @@ namespace SEB_NAMESPACE {
     }
   }
   
-  template<typename Float>
-  void Subspan<Float>::hessenberg_clear (unsigned int pos)
+  template<typename Float, class Pt, class PointAccessor>
+  void Subspan<Float, Pt, PointAccessor>::hessenberg_clear (unsigned int pos)
   // Given R in lower Hessenberg form with subdiagonal entries 0 to
   // pos-1 already all zero, clears the remaining subdiagonal entries
   // via Givens rotations.
@@ -309,8 +311,8 @@ namespace SEB_NAMESPACE {
     }
   }
   
-  template<typename Float>
-  void Subspan<Float>::special_rank_1_update ()
+  template<typename Float, class Pt, class PointAccessor>
+  void Subspan<Float, Pt, PointAccessor>::special_rank_1_update ()
   // Update current QR-decomposition "A = QR" to
   // A + u * [1,...,1] = Q' R'.
   {
diff --git a/cpp/main/Subspan.h b/cpp/main/Subspan.h
index a18584a..0d8f62b 100644
--- a/cpp/main/Subspan.h
+++ b/cpp/main/Subspan.h
@@ -6,8 +6,6 @@
 #ifndef SEB_SUBSPAN_H
 #define SEB_SUBSPAN_H
 
-#include "Seb_point.h"
-
 namespace SEB_NAMESPACE {
   
   template<typename Float>
@@ -16,7 +14,7 @@ namespace SEB_NAMESPACE {
     return x * x;
   }
   
-  template<typename Float>
+  template<typename Float, class Pt, class PointAccessor>
   class Subspan
   // An instance of this class represents the affine hull of a
   // non-empty set M of affinely independent points.  The set M is not
@@ -73,12 +71,9 @@ namespace SEB_NAMESPACE {
   //   in other words, to the point in S with index global_index(i).
   //   Complexity: O(dim^2)
   {
-  private: // types:
-    typedef Point<Float> Pt;
-    
   public: // construction and deletion:
     
-    Subspan(unsigned int dim, const std::vector<Pt>& S, int i);
+    Subspan(unsigned int dim, const PointAccessor& S, int i);
     // Constructs an instance representing the affine hull aff(M) of M={p},
     // where p is the point S[i] from S.
     //
@@ -154,7 +149,7 @@ namespace SEB_NAMESPACE {
     // A + u [1,...,1] = Q' R'.
     
   private: // member fields:
-    const std::vector<Pt>& S;          // a const-reference to the set S
+    const PointAccessor &S;            // a const-reference to the set S
     std::vector<bool> membership;      // S[i] in M iff membership[i]
     const unsigned int dim;            // ambient dimension (not to be
     // confused with the rank r,
diff --git a/cpp/test/example.C b/cpp/test/example.C
index f9679a4..6ac40e7 100644
--- a/cpp/test/example.C
+++ b/cpp/test/example.C
@@ -6,12 +6,13 @@
 #include <iostream>
 #include <cstdio>
 
-#include <Seb.h>
-#include <Seb_debug.C> // ... only because we use Seb::Timer below
+#include "Seb.h"
+#include "Seb_debug.C" // ... only because we use Seb::Timer below
 
 int main(int argn,char **argv) {
   typedef double FT;
   typedef Seb::Point<FT> Point;
+  typedef std::vector<Point> PointVector;
   typedef Seb::Smallest_enclosing_ball<FT> Miniball;
 
   using std::cout;
@@ -31,7 +32,7 @@ int main(int argn,char **argv) {
   const int n = std::atoi(argv[1]), d = std::atoi(argv[2]);
 
   // construct n random points in dimension d:
-  vector<Point> S;
+  PointVector S;
   vector<double> coords(d);
   srand(clock());
   for (int i=0; i<n; ++i) {
@@ -55,10 +56,8 @@ int main(int argn,char **argv) {
        << "====================================================" << endl;
   Seb::Timer::instance().start("all");
 
-  // compute the miniball:
-  Miniball mb(d);
-  for (unsigned int i=0; i<S.size(); ++i)
-    mb.insert(S[i].begin());
+  // compute the miniball by inserting each value
+  Miniball mb(d, S);
 
   // output:
   FT rad = mb.radius();