docs: public Clustering<SubspaceModel> run(Database database, Relation<V> relation) {
      if(RelationUtil.dimensionality(relation) < l) {
        throw new IllegalStateException("Dimensionality of data < parameter l! (" + RelationUtil.dimensionality(relation) + " < " + l + ")");
      }
      DistanceQuery<V> distFunc = database.getDistanceQuery(relation, SquaredEuclideanDistanceFunction.STATIC);
      RangeQuery<V> rangeQuery = database.getRangeQuery(distFunc);
      final Random random = rnd.getSingleThreadedRandom();
  
      // TODO: use a StepProgress!
      // initialization phase
      if(LOG.isVerbose()) {
        LOG.verbose("1. Initialization phase...");
      }
      int sampleSize = Math.min(relation.size(), k_i * k);
      DBIDs sampleSet = DBIDUtil.randomSample(relation.getDBIDs(), sampleSize, random);
  
      int medoidSize = Math.min(relation.size(), m_i * k);
      ArrayDBIDs medoids = greedy(distFunc, sampleSet, medoidSize, random);
  
      if(LOG.isDebugging()) {
        LOG.debugFine(new StringBuilder().append("sampleSize ").append(sampleSize).append('\n') //
            .append("sampleSet ").append(sampleSet).append('\n') //
            .append("medoidSize ").append(medoidSize).append('\n') //
            .append("m ").append(medoids).toString());
      }
  
      // iterative phase
      if(LOG.isVerbose()) {
        LOG.verbose("2. Iterative phase...");
      }
      double bestObjective = Double.POSITIVE_INFINITY;
      ArrayDBIDs m_best = null;
      DBIDs m_bad = null;
      ArrayDBIDs m_current = initialSet(medoids, k, random);
  
      if(LOG.isDebugging()) {
        LOG.debugFine(new StringBuilder().append("m_c ").append(m_current).toString());
      }
  
      IndefiniteProgress cprogress = LOG.isVerbose() ? new IndefiniteProgress("Current number of clusters:", LOG) : null;
  
      ArrayList<PROCLUSCluster> clusters = null;
      int loops = 0;
      while(loops < 10) {
        long[][] dimensions = findDimensions(m_current, relation, distFunc, rangeQuery);
        clusters = assignPoints(m_current, dimensions, relation);
        double objectiveFunction = evaluateClusters(clusters, dimensions, relation);
  
        if(objectiveFunction < bestObjective) {
          // restart counting loops
          loops = 0;
          bestObjective = objectiveFunction;
          m_best = m_current;
          m_bad = computeBadMedoids(m_current, clusters, (int) (relation.size() * 0.1 / k));
        }
  
        m_current = computeM_current(medoids, m_best, m_bad, random);
        loops++;
        if(cprogress != null) {
          cprogress.setProcessed(clusters.size(), LOG);
        }
      }
  
      LOG.setCompleted(cprogress);
  
      // refinement phase
      if(LOG.isVerbose()) {
        LOG.verbose("3. Refinement phase...");
      }
  
      List<Pair<double[], long[]>> dimensions = findDimensions(clusters, relation);
      List<PROCLUSCluster> finalClusters = finalAssignment(dimensions, relation);
  
      // build result
      int numClusters = 1;
      Clustering<SubspaceModel> result = new Clustering<>("ProClus clustering", "proclus-clustering");
      for(PROCLUSCluster c : finalClusters) {
        Cluster<SubspaceModel> cluster = new Cluster<>(c.objectIDs);
        cluster.setModel(new SubspaceModel(new Subspace(c.getDimensions()), c.centroid));
        cluster.setName("cluster_" + numClusters++);
  
        result.addToplevelCluster(cluster);
      }
      return result;
    }