best_partition.cpp

#include <Rcpp.h>
using namespace Rcpp;

double distance(double x, double y) {
  double temp = x - y;
  return temp * temp;
}

// [[Rcpp::export]]
int bestLeftEdge(NumericVector sorted_values) {
  int value_count = sorted_values.size();

  int left_edge = value_count / 2;
  int left_count = left_edge + 1;

  int right_edge = left_edge + 1;
  int right_count = value_count - left_count;

  double left_total = 0.0;
  for (int i = 0; i <= left_edge; i++) {
    left_total += sorted_values[i];
  }
  double left_mean = left_total / left_count;

  double right_total = 0.0;
  for (int i = right_edge; i < value_count; i++) {
    right_total += sorted_values[i];
  }
  double right_mean = right_total / right_count;

  int loop = 1;
  double left_value;
  double right_value;

  while(loop == 1) {
    left_value = sorted_values[left_edge];
    right_value = sorted_values[right_edge];

    if (distance(left_value, right_mean) < distance(left_value, left_mean)) {
      right_total += left_value;
      right_count++;
      right_mean = right_total / right_count;
      right_edge--;

      left_total -= left_value;
      left_count--;
      left_mean = left_total / left_count;
      left_edge--;
    } else if (distance(right_value, left_mean) < distance(right_value, right_mean)) {
      right_total -= right_value;
      right_count--;
      right_mean = right_total / right_count;
      right_edge++;

      left_total += right_value;
      left_count++;
      left_mean = left_total / left_count;
      left_edge++;
    } else {
      loop = 0;
    }   
  }

  return left_edge + 1;
}