Initial commit.

.gitignore

+*.netbeans*
+*.xml
+*.properties
+*.class
+*.mf

src/solvers/Greedy.java

+package solvers;
+
+import java.util.ArrayList;
+import java.util.Collection;
+import java.util.Collections;
+import java.util.HashSet;
+import java.util.List;
+import xmtrclusters.Problem;
+
+/**
+ * Greedy implements a greedy heuristic for the clustering problem.
+ *
+ * Initially, each transmitter is a cluster unto itself. In each iteration,
+ * all legal merges of a pair of clusters are tested, and the one with maximum
+ * benefit is implemented. If no merges are beneficial and the number of
+ * clusters exceeds the maximum allowed, the least harmful legal merger is done.
+ * This repeats until the solution is legal and cannot be improved by a merger,
+ * or the solution violates the cluster limit and no pair of clusters can
+ * be merged.
+ *
+ * @author Paul A. Rubin (rubin@msu.edu)
+ */
+
+record Merger(int first, int second, double change) { }
+
+public final class Greedy {
+  private final List<Collection<Integer>> clusters;  // list of clusters
+  private final Problem problem;                     // problem to solve
+  private final int maxSize;                         // maximum problem size
+
+
+  /**
+   * Constructs the heuristic instance and solves it.
+   * @param prob the problem to solve
+   */
+  public Greedy(final Problem prob) {
+    long time = System.currentTimeMillis();
+    problem = prob;
+    maxSize = prob.getMaxSize();           // max cluster size
+    int maxCount = prob.getMaxClusters();  // max cluster count
+    // Initialize the cluster list with singleton clusters.
+    clusters = new ArrayList<>();
+    for (int t = 0; t < prob.getNTrans(); t++) {
+      HashSet<Integer> c = new HashSet<>();
+      c.add(t);
+      clusters.add(c);
+    }
+    // Loop indefinitely until no further changes can be made.
+    while (true) {
+      // Find the pair of clusters that can be merged with maximum benefit.
+      Merger merger = findMerger();
+      // If no merger is possible due to size limits, end the loop.
+      if (merger == null) {
+        break;
+      }
+      // If a beneficial merger was found, make it and loop.
+      if (merger.change() > 0) {
+        execute(merger);
+        continue;
+      }
+      // If we fall to here, no beneficial merge is possible. Check compliance
+      // with the cluster count limit.
+      if (clusters.size() <= maxCount) {
+        // The solution is in compliance, so exit.
+        break;
+      } else {
+        // If a non-beneficial merger was found, do it and continue looping.
+        if (merger == null) {
+          // Out of compliance with no possible merger.
+          clusters.clear();
+          break;
+        } else {
+          execute(merger);
+        }
+      }
+    }
+    time = System.currentTimeMillis() - time;
+    System.out.println("\nThe greedy heuristic ran for " + time + " ms.");
+  }
+
+  /**
+   * Finds the best feasible merger given the current list of clusters.
+   * @return a merger record showing the indices of two clusters to merge and
+   * the change in objective value, or null if there are no legal mergers
+   */
+  private Merger findMerger() {
+    int first = -1;
+    int second = -1;
+    double benefit = Double.MIN_VALUE;
+    for (int i = 0; i < clusters.size(); i++) {
+      Collection<Integer> c1 = clusters.get(i);
+      double obj1 = problem.clusterQuality(c1);
+      for (int j = i + 1; j < clusters.size(); j++) {
+        Collection<Integer> c2 = clusters.get(j);
+        // Make sure the size limit is respected.
+        if (c1.size() + c2.size() <= maxSize) {
+          // Compute the benefit of a merge.
+          Collection<Integer> z = new HashSet<>(c1);
+          z.addAll(c2);
+          double delta =
+            problem.clusterQuality(z) - obj1 - problem.clusterQuality(c2);
+          if (delta > benefit) {
+            benefit = delta;
+            first = i;
+            second = j;
+          }
+        }
+      }
+    }
+    // Return the merger (or null if none was found).
+    if (first < 0) {
+      return null;
+    } else {
+      return new Merger(first, second, benefit);
+    }
+  }
+
+  /**
+   * Implements a merger.
+   * @param merger the merger to implement
+   */
+  private void execute(final Merger merger) {
+    HashSet<Integer> merged = new HashSet<>(clusters.get(merger.first()));
+    merged.addAll(clusters.get(merger.second()));
+    clusters.remove(merger.second());
+    clusters.remove(merger.first());
+    clusters.add(merged);
+  }
+
+  /**
+   * Gets the final result of the heuristic.
+   * @return the collection of clusters found
+   */
+  public Collection<Collection<Integer>> getSolution() {
+    return Collections.unmodifiableCollection(clusters);
+  }
+}

src/solvers/MIPModel.java

+package solvers;
+
+import ilog.concert.IloException;
+import ilog.concert.IloIntVar;
+import ilog.concert.IloLinearNumExpr;
+import ilog.concert.IloNumVar;
+import ilog.cplex.IloCplex;
+import java.util.ArrayList;
+import java.util.Collection;
+import java.util.HashSet;
+import xmtrclusters.Problem;
+
+/**
+ * MIPModel implements a mixed-integer linear program for the clustering
+ * problem.
+ *
+ * @author Paul A. Rubin (rubin@msu.edu)
+ */
+public final class MIPModel {
+  private static final double HALF = 0.5;  // for rounding binary values
+
+  private final Problem problem;  // the parent problem
+
+  // Model components.
+  private final IloCplex mip;     // the MIP model
+  private final IloIntVar[][] x;  // x[i][j] = 1 if trans i & j are in the same
+                                  // cluster; x[i][i] = 1 if i anchors a cluster
+  private final IloIntVar[][] y;  // y[i][j] = 1 if trans i & user j are
+                                  // in the same cluster
+  private final IloNumVar[][] z;  // z[i][j] is the objective contribution of
+                                  // the trans i / user j pairing
+  private final IloNumVar[] q;    // q[j] is the solution quality for user j
+  private final IloNumVar[][] v;  // v[i][j] = 1 if transmitters i and j
+                                  // are in the same cluster
+
+  /**
+   * Constructor.
+   * @param prob the problem to solve
+   * @throws IloException if CPLEX throws an error building the model
+   */
+  public MIPModel(final Problem prob) throws IloException {
+    problem = prob;
+    int nU = problem.getNUsers();
+    int nT = problem.getNTrans();
+    double[] qMax = problem.getMaxQuality();
+    // Initialize the model.
+    mip = new IloCplex();
+    // Add the variables.
+    x = new IloIntVar[nT][nT];
+    y = new IloIntVar[nT][nU];
+    z = new IloNumVar[nT][nU];
+    q = new IloNumVar[nU];
+    v = new IloNumVar[nT][nT];
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = 0; t1 < nT; t1++) {
+        x[t][t1] = mip.boolVar("x_" + t + "_" + t1);
+        v[t][t1] = mip.numVar(0, 1, "v_" + t + "_" + t1);
+      }
+      for (int u = 0; u < nU; u++) {
+        y[t][u] = mip.boolVar("y_" + t + "_" + u);
+        z[t][u] = mip.numVar(0, qMax[u], "z_" + t + "_" + u);
+      }
+    }
+    for (int u = 0; u < nU; u++) {
+      q[u] = mip.numVar(0, qMax[u], "q_" + u);
+    }
+    // The objective is to maximize the total userQuality.
+    mip.addMaximize(mip.sum(q));
+    // Make the x variables symmetric.
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = t + 1; t1 < nT; t1++) {
+        mip.addEq(x[t][t1], x[t1][t], "symmetry_" + t + "_" + t1);
+      }
+    }
+    // Limit the number of anchors (hence the number of clusters).
+    IloLinearNumExpr expr = mip.linearNumExpr();
+    for (int t = 0; t < nT; t++) {
+      expr.addTerm(1.0, x[t][t]);
+    }
+    mip.addLe(expr, problem.getMaxClusters(), "max_cluster_count");
+    // Limit the size of every cluster.
+    int maxS = problem.getMaxSize();
+    for (int t = 0; t < nT; t++) {
+      mip.addLe(mip.diff(mip.sum(x[t]), x[t][t]), maxS - 1,
+                "cluster_size_" + t);
+    }
+    // Avoid multiple anchors in a cluster.
+    // (These constraints may be unnecessary.)
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = t + 1; t1 < nT; t1++) {
+        mip.addLe(mip.sum(x[t][t], x[t][t1], x[t1][t1]), 2.0,
+                  "one_anchor_" + t + "_" + t1);
+      }
+    }
+    // Force each transmitter to be in a cluster with an anchor.
+    for (int t = 0; t < nT; t++) {
+      expr = mip.linearNumExpr();
+      expr.addTerm(1.0, x[t][t]);
+      for (int t1 = 0; t1 < nT; t1++) {
+        if (t != t1) {
+          mip.addLe(v[t][t1], x[t][t1], "v_def1_" + t + "_" + t1);
+          mip.addLe(v[t][t1], x[t1][t1], "v_def2_" + t + "_" + t1);
+          expr.addTerm(1.0, v[t][t1]);
+        }
+      }
+      mip.addEq(expr, 1.0, "in_cluster_" + t);
+    }
+    // Cluster membership is transitive.
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = 0; t1 < nT; t1++) {
+        if (t != t1) {
+          for (int t2 = 0; t2 < nT; t2++) {
+            if (t2 != t && t2 != t1) {
+              mip.addGe(x[t][t2], mip.diff(mip.sum(x[t][t1], x[t1][t2]), 1.0),
+                        "transitive_" + t + "_" + t1 + "_" + t2);
+            }
+          }
+        }
+      }
+    }
+    // Antisymmetry: lowest index transmitter is the anchor.
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = t + 1; t1 < nT; t1++) {
+        mip.addLe(x[t1][t1], mip.diff(1.0, x[t][t1]),
+                  "asymmetry_" + t + "_" + t1);
+      }
+    }
+    // Put each user in with their favorite transmitter.
+    int[] tau = problem.getBest();
+    for (int u = 0; u < nU; u++) {
+      mip.addEq(y[tau[u]][u], 1.0, "favorite_transmitter_" + u);
+    }
+    // Determine other transmitter-user pairings.
+    for (int u = 0; u < nU; u++) {
+      for (int t = 0; t < nT; t++) {
+        if (t != tau[u]) {
+          mip.addEq(y[t][u], x[tau[u]][t], "pairing_" + t + "_" + u);
+        }
+      }
+    }
+    // Define the z variables.
+    for (int t = 0; t < nT; t++) {
+      for (int u = 0; u < nU; u++) {
+        mip.addLe(z[t][u], mip.prod(qMax[u], y[t][u]), "zdef1_" + t + "_" + u);
+        expr = mip.linearNumExpr();
+        expr.addTerm(1.0, z[t][u]);
+        expr.addTerm(-1.0, q[u]);
+        expr.addTerm(qMax[u], y[t][u]);
+        mip.addLe(expr, qMax[u], "zdef2_" + t + "_" + u);
+        expr = mip.linearNumExpr();
+        expr.addTerm(1.0, z[t][u]);
+        expr.addTerm(-1.0, q[u]);
+        expr.addTerm(-qMax[u], y[t][u]);
+        mip.addGe(expr, -qMax[u], "zdef3_" + t + "_" + u);
+      }
+    }
+    // Linearize q.
+    for (int u = 0; u < nU; u++) {
+      expr = mip.linearNumExpr();
+      double[] w = problem.getWeights(u);
+      for (int i = 0; i < nT; i++) {
+        expr.addTerm(w[i], q[u]);
+        expr.addTerm(-w[i], z[i][u]);
+        expr.addTerm(-w[i], y[i][u]);
+      }
+      mip.addEq(expr, 0.0, "qdef_" + u);
+    }
+  }
+
+  /**
+   * Exports the model to a file.
+   * @param target the target file path/name
+   * @throws IloException if CPLEX encounters an error during export
+   */
+  public void export(final String target) throws IloException {
+    mip.exportModel(target);
+  }
+
+  /**
+   * Solves the model.
+   * @param sec time limit in seconds
+   * @return the final objective value
+   * @throws IloException if CPLEX encounters an error
+   */
+  public double solve(final double sec) throws IloException {
+    mip.setParam(IloCplex.DoubleParam.TimeLimit, sec);
+    mip.setParam(IloCplex.Param.Emphasis.MIP, 3);
+    mip.solve();
+    return mip.getObjValue();
+  }
+
+  /**
+   * Gets the transmitter clusters in the final solution.
+   * @return the collection of clusters
+   * @throws IloException if CPLEX cannot extract the solution.
+   */
+  public Collection<Collection<Integer>> getClusters() throws IloException {
+    // Make sure a solution exists.
+    IloCplex.Status status = mip.getStatus();
+    if (status != IloCplex.Status.Optimal
+        && status != IloCplex.Status.Feasible) {
+      throw new IllegalArgumentException("Cannot extract a"
+                                         + " nonexistent solution!");
+    }
+    int nT = problem.getNTrans();
+    // Extract the values of x[][] in the solution.
+    boolean[][] xx = new boolean[nT][nT];
+    for (int t = 0; t < nT; t++) {
+      for (int t1 = t; t1 < nT; t1++) {
+        xx[t][t1] = mip.getValue(x[t][t1]) > HALF;
+        xx[t1][t] = xx[t][t1];
+      }
+    }
+    // Identify the clusters.
+    ArrayList<Collection<Integer>> clusters = new ArrayList<>();
+    for (int t = 0; t < nT; t++) {
+      // Look for anchors.
+      if (xx[t][t]) {
+        // t is an anchor; add the cluster containing it (ncluding itself).
+        HashSet<Integer> c = new HashSet<>();
+        for (int t1 = 0; t1 < nT; t1++) {
+          if (xx[t][t1]) {
+            c.add(t1);
+          }
+        }
+        clusters.add(c);
+      }
+    }
+    return clusters;
+  }
+}

src/solvers/MIPModel2.java

+package solvers;
+
+import ilog.concert.IloException;
+import ilog.concert.IloIntVar;
+import ilog.concert.IloLinearNumExpr;
+import ilog.concert.IloNumVar;
+import ilog.cplex.IloCplex;
+import java.util.ArrayList;
+import java.util.Collection;
+import java.util.HashSet;
+import xmtrclusters.Problem;
+
+/**
+ * MIPModel2 implements a mixed-integer programming model for the clustering
+ * problem.
+ *
+ * It uses the same linearization of the objective that MIPModel uses, but
+ * handles cluster creation as an assignment model with some symmetry-breaking
+ * constraints.
+ *
+ * The model winds up being considerably larger than MIPModel and takes several
+ * times as long to solve.
+ *
+ * @author Paul A. Rubin (rubin@msu.edu)
+ */
+public final class MIPModel2 {
+  private static final double HALF = 0.5;  // for rounding binary values
+
+  private final Problem problem;  // the parent problem
+
+  // Model components.
+  private final IloCplex mip;     // the MIP model
+  private final IloIntVar[][] x;  // x[i][j] = 1 if trans i & j are in the same
+                                  // cluster; x[i][i] = 1 if i anchors a cluster
+  private final IloIntVar[][] y;  // y[i][j] = 1 if trans i & user j are
+                                  // in the same cluster
+  private final IloNumVar[][] z;  // z[i][j] is the objective contribution of
+                                  // the trans i / user j pairing
+  private final IloNumVar[] q;    // q[j] is the solution quality for user j
+  private final IloNumVar[][][] v;
+    // v[t][u][c] is an indicator for transmitter t and user u both belonging
+    // to cluster c
+
+  /**
+   * Constructor.
+   * @param prob the problem to solve
+   * @throws IloException if CPLEX throws an error building the model
+   */
+  public MIPModel2(final Problem prob) throws IloException {
+    problem = prob;
+    int nU = problem.getNUsers();
+    int nT = problem.getNTrans();
+    int nC = problem.getMaxClusters();
+    int maxC = problem.getMaxSize();
+    double[] qMax = problem.getMaxQuality();
+    // Initialize the model.
+    mip = new IloCplex();
+    // Add the variables.
+    x = new IloIntVar[nT][nC];
+    y = new IloIntVar[nT][nU];
+    z = new IloNumVar[nT][nU];
+    q = new IloNumVar[nU];
+    v = new IloNumVar[nT][nU][nC];
+    for (int t = 0; t < nT; t++) {
+      for (int c = 0; c < nC; c++) {
+        x[t][c] = mip.boolVar("x_" + t + "_" + c);
+      }
+      for (int u = 0; u < nU; u++) {
+        y[t][u] = mip.boolVar("y_" + t + "_" + u);
+        z[t][u] = mip.numVar(0, qMax[u], "z_" + t + "_" + u);
+        for (int c = 0; c < nC; c++) {
+          v[t][u][c] = mip.numVar(0, 1, "v_" + t + "_" + u + "_" + c);
+        }
+      }
+    }
+    for (int u = 0; u < nU; u++) {
+      q[u] = mip.numVar(0, qMax[u], "q_" + u);
+    }
+    // The objective is to maximize the total userQuality.
+    mip.addMaximize(mip.sum(q));
+    // Every transmitter is assigned to exactly one cluster.
+    for (int t = 0; t < nT; t++) {
+      mip.addEq(mip.sum(x[t]), 1.0, "assign_" + t + "_once");
+    }
+    // Observe cluster capacity limits.
+    for (int c = 0; c < nC; c++) {
+      IloLinearNumExpr expr = mip.linearNumExpr();
+      for (int t = 0; t < nT; t++) {
+        expr.addTerm(1.0, x[t][c]);
+      }
+      mip.addLe(expr, maxC, "capacity_" + c);
+    }
+    // Antisymmetry: Automatically assign transmitter 0 to cluster 0.
+    x[0][0].setLB(1.0);
+    // Antisymmetry: To assign transmitter t to cluster c, cluster c-1 must
+    // contain at least one transmitter with index less than t.
+    for (int t = 1; t < nT; t++) {
+      for (int c = 1; c < nC; c++) {
+        IloLinearNumExpr expr = mip.linearNumExpr();
+        for (int t0 = 0; t0 < t; t0++) {
+          expr.addTerm(1.0, x[t0][c - 1]);
+        }
+        mip.addLe(x[t][c], expr, "asym_" + t + "_" + c);
+      }
+    }
+    // Every user is assigned to the cluster containing its best transmitter.
+    int[] tau = problem.getBest();
+    for (int u = 0; u < nU; u++) {
+      mip.addEq(y[tau[u]][u], 1.0, "favorite_transmitter_" + u);
+    }
+    // Use the auxiliary (v) variables to define y[t][u] when t is not tau[u].
+    for (int t = 0; t < nT; t++) {
+      for (int u = 0; u < nU; u++) {
+        IloLinearNumExpr expr = mip.linearNumExpr();
+        for (int c = 0; c < nC; c++) {
+          // v[t][u][c] = 1 if t is in c and if u is in c (because tau[u] is
+          // in c).
+          mip.addLe(v[t][u][c], x[tau[u]][c], "v1_" + t + "_" + u + "_" + c);
+          mip.addLe(v[t][u][c], x[t][c], "v2_" + t + "_" + u + "_" + c);
+          mip.addGe(v[t][u][c], mip.diff(mip.sum(x[tau[u]][c], x[t][c]), 1.0),
+                    "v3_" + t + "_" + u + "_" + c);
+          expr.addTerm(1.0, v[t][u][c]);
+        }
+        // y[t][u] = 1 if v[t][u][c] = 1 for some c.
+        mip.addEq(y[t][u], expr, "def_y_" + t + "_" + u);
+      }
+    }
+    // Define the z variables.
+    for (int t = 0; t < nT; t++) {
+      for (int u = 0; u < nU; u++) {
+        mip.addLe(z[t][u], mip.prod(qMax[u], y[t][u]), "zdef1_" + t + "_" + u);
+        IloLinearNumExpr expr = mip.linearNumExpr();
+        expr.addTerm(1.0, z[t][u]);
+        expr.addTerm(-1.0, q[u]);
+        expr.addTerm(qMax[u], y[t][u]);
+        mip.addLe(expr, qMax[u], "zdef2_" + t + "_" + u);
+        expr = mip.linearNumExpr();
+        expr.addTerm(1.0, z[t][u]);
+        expr.addTerm(-1.0, q[u]);
+        expr.addTerm(-qMax[u], y[t][u]);
+        mip.addGe(expr, -qMax[u], "zdef3_" + t + "_" + u);
+      }
+    }
+    // Linearize q.
+    for (int u = 0; u < nU; u++) {
+      IloLinearNumExpr expr = mip.linearNumExpr();
+      double[] w = problem.getWeights(u);
+      for (int i = 0; i < nT; i++) {
+        expr.addTerm(w[i], q[u]);
+        expr.addTerm(-w[i], z[i][u]);
+        expr.addTerm(-w[i], y[i][u]);
+      }
+      mip.addEq(expr, 0.0, "qdef_" + u);
+    }
+  }
+
+  /**
+   * Exports the model to a file.
+   * @param target the target file path/name
+   * @throws IloException if CPLEX encounters an error during export
+   */
+  public void export(final String target) throws IloException {
+    mip.exportModel(target);
+  }
+
+  /**
+   * Solves the model.
+   * @param sec time limit in seconds
+   * @return the final objective value
+   * @throws IloException if CPLEX encounters an error
+   */
+  public double solve(final double sec) throws IloException {
+    mip.setParam(IloCplex.DoubleParam.TimeLimit, sec);
+    mip.setParam(IloCplex.Param.Emphasis.MIP, 2);
+    mip.solve();
+    return mip.getObjValue();
+  }
+  /**
+   * Gets the transmitter clusters in the final solution.
+   * @return the collection of clusters
+   * @throws IloException if CPLEX cannot extract the solution.
+   */
+  public Collection<Collection<Integer>> getClusters() throws IloException {
+    // Make sure a solution exists.
+    IloCplex.Status status = mip.getStatus();
+    if (status != IloCplex.Status.Optimal
+        && status != IloCplex.Status.Feasible) {
+      throw new IllegalArgumentException("Cannot extract a"
+                                         + " nonexistent solution!");
+    }
+    int nT = problem.getNTrans();
+    int nC = problem.getMaxClusters();
+    // Extract the values of x[][] in the solution.
+    boolean[][] xx = new boolean[nT][nC];
+    for (int t = 0; t < nT; t++) {
+      for (int c = 0; c < nC; c++) {
+        xx[t][c] = mip.getValue(x[t][c]) > HALF;
+      }
+    }
+    // Identify the clusters.
+    ArrayList<Collection<Integer>> clusters = new ArrayList<>();
+    for (int c = 0; c < nC; c++) {
+      HashSet<Integer> cluster = new HashSet<>();
+      for (int t = 0; t < nT; t++) {
+        if (xx[t][c]) {
+          cluster.add(t);
+        }
+      }
+      clusters.add(cluster);
+    }
+    return clusters;
+  }
+}

src/xmtrclusters/Problem.java

+package xmtrclusters;
+
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.Collection;
+import java.util.HashMap;
+import java.util.HashSet;
+import java.util.List;
+import java.util.Random;
+import java.util.Set;
+import java.util.stream.Collectors;
+import java.util.stream.IntStream;
+
+/**
+ * Problem holds the data for a problem instance.
+ *
+ * @author Paul A. Rubin (rubin@msu.edu)
+ */
+public final class Problem {
+  private final int nTrans;         // # of transmitters
+  private final int nUsers;         // # of users
+  private final int maxClusters;    // maximum # of clusters
+  private final int maxSize;        // maximum # of transmitters per cluster
+  private final double[][] weight;  // weight for transmitter-user pair
+  private final int[] best;         // index of best transmitter for each user
+  private final List<Integer> allTransmitters;  // list of all transmitters
+  private final double[] maxQuality;  // maximum possible user for each user
+  private final HashMap<Integer, Set<Integer>> assignedUsers;
+    // maps each transmitter to the users best served by it
+
+  /**
+   * Constructs a random problem instance.
+   * @param nT the number of transmitters
+   * @param nU the number of users
+   * @param maxC the maximum number of clusters
+   * @param maxT the maximum number of transmitters
+   * @param seed the seed for the random number generator.
+   */
+  public Problem(final int nT, final int nU, final int maxC, final int maxT,
+                 final int seed) {
+    // Set the problem dimensions.
+    nTrans = nT;
+    nUsers = nU;
+    maxClusters = maxC;
+    maxSize = maxT;
+    weight = new double[nTrans][nUsers];
+    best = new int[nUsers];
+    allTransmitters = IntStream.range(0, nTrans).boxed()
+                               .collect(Collectors.toList());
+    assignedUsers = new HashMap<>();
+    for (int t = 0; t < nTrans; t++) {
+      assignedUsers.put(t, new HashSet<>());
+    }
+    // Assign service quality values (weights) randomly and determine the
+    // best transmitter for each user.
+    Random rng = new Random(seed);
+    for (int u = 0; u < nUsers; u++) {
+      double z = -1;
+      for (int t = 0; t < nTrans; t++) {
+        weight[t][u] = rng.nextDouble();
+        if (weight[t][u] > z) {
+          best[u] = t;
+          z = weight[t][u];
+        }
+      }
+      assignedUsers.get(best[u]).add(u);
+    }
+    // The maximum quality for each user occurs when the worst transmitter for
+    // that user is the only transmitter not in the user's cluster.
+    maxQuality = new double[nUsers];
+    for (int u = 0; u < nUsers; u++) {
+      double worst = Double.MAX_VALUE;
+      int w = -1;
+      for (int t = 0; t < nTrans; t++) {
+        if (weight[t][u] < worst) {
+          worst = weight[t][u];
+          w = t;
+        }
+      }
+      HashSet<Integer> c = new HashSet<>(allTransmitters);
+      c.remove(w);
+      maxQuality[u] = userQuality(u, c);
+    }
+  }
+
+  /**
+   * Computes the service quality for a single user.
+   * @param user the index of the user
+   * @param cluster the indices of all transmitters in the same cluster as the
+   * user
+   * @return the service quality for the user
+   */
+  public double userQuality(final int user, final Collection<Integer> cluster) {
+    // Make sure at least one transmitter is not in the cluster (to avoid
+    // division by zero).
+    if (cluster.size() == nTrans) {
+      throw new IllegalArgumentException("Cannot have all transmitters"
+                                         + " in one cluster.");
+    }
+    // The numerator is the sum of weights for all transmitters in the cluster.
+    double num = cluster.stream().mapToDouble(t -> weight[t][user]).sum();
+    // The denominator is the sum of weights for all transmitters not in
+    // the cluster.
+    HashSet<Integer> out = new HashSet<>(allTransmitters);
+    out.removeAll(cluster);
+    double denom = out.stream().mapToDouble(t -> weight[t][user]).sum();
+    return num / denom;
+  }
+
+  /**
+   * Generates a summary of a solution.
+   * @param solution the solution (a collection of clusters of transmitters)
+   * @return a string summarizing the solution.
+   */
+  public String report(final Collection<Collection<Integer>> solution) {
+    StringBuilder sb = new StringBuilder();
+    sb.append("The solution has ").append(solution.size())
+      .append(" clusters.\n");
+    for (Collection<Integer> cluster : solution) {
+      sb.append("\tCluster ")
+        .append(Arrays.toString(cluster.stream().sorted()
+                                       .mapToInt(i -> i).toArray()))
+        .append(" serves ");
+      ArrayList<Integer> list = new ArrayList<>();
+      cluster.stream().forEach(t -> list.addAll(assignedUsers.get(t)));
+      sb.append(Arrays.toString(list.stream().sorted()
+                                    .mapToInt(i -> i).toArray()))
+        .append("\n");
+    }
+    sb.append("\nOverall quality = ").append(totalQuality(solution))
+                                     .append(".");
+    return sb.toString();
+  }
+
+  /**
+   * Gets the number of transmitters.
+   * @return the number of transmitters
+   */
+  public int getNTrans() {
+    return nTrans;
+  }
+
+  /**
+   * Gets the number of users.
+   * @return the number of users
+   */
+  public int getNUsers() {
+    return nUsers;
+  }
+
+  /**
+   * Gets the maximum number of clusters.
+   * @return the maximum number of clusters
+   */
+  public int getMaxClusters() {
+    return maxClusters;
+  }
+
+  /**
+   * Gets the maximum cluster size.
+   * @return the maximum cluster size
+   */
+  public int getMaxSize() {
+    return maxSize;
+  }
+
+  /**
+   * Gets the maximum possible quality for each user.
+   * @return the vector of quality bounds
+   */
+  public double[] getMaxQuality() {
+    return Arrays.copyOf(maxQuality, maxQuality.length);
+  }
+
+  /**
+   * Gets the indices of the best transmitters for each user.
+   * @return the vector of best transmitter indices
+   */
+  public int[] getBest() {
+    return Arrays.copyOf(best, best.length);
+  }
+
+  /**
+   * Gets the weights for a given user.
+   * @param user the index of the user
+   * @return the vector of weights for all transmitters and that user
+   */
+  public double[] getWeights(final int user) {
+    return allTransmitters.stream().mapToDouble(t -> weight[t][user]).toArray();
+  }
+
+  /**
+   * Computes the total quality of a proposed solution.
+   * @param clusters a collection of clusters
+   * @return the total user quality for that solution
+   */
+  public double totalQuality(final Collection<Collection<Integer>> clusters) {
+    return clusters.stream().mapToDouble(c -> clusterQuality(c)).sum();
+  }
+
+  /**
+   * Computes the total quality measure for all users in a cluster.
+   * @param cluster the cluster of transmitters
+   * @return the total quality for users in that cluster
+   */
+  public double clusterQuality(final Collection<Integer> cluster) {
+    // Find the users in the cluster.
+    ArrayList<Integer> list = new ArrayList<>();
+    cluster.stream().forEach(t -> list.addAll(assignedUsers.get(t)));
+    // Total their quality values.
+    return list.stream().mapToDouble(u -> userQuality(u, cluster)).sum();
+  }
+}

src/xmtrclusters/XmtrClusters.java

+package xmtrclusters;
+
+import ilog.concert.IloException;
+import java.util.Collection;
+import solvers.Greedy;
+import solvers.MIPModel;
+import solvers.MIPModel2;
+
+/**
+ * XmtrClusters experiments with solutions to a problem of clustering
+ * transmitters and users.
+ *
+ * Problem source: https://or.stackexchange.com/questions/7471/how-to-perform-
+ * clustering-of-two-different-sets-of-entities/
+ *
+ * @author Paul A. Rubin (rubin@msu.edu)
+ */
+public final class XmtrClusters {
+
+  /**
+   * Dummy constructor.
+   */
+  private XmtrClusters() { }
+
+  /**
+   * Runs the experiments.
+   * @param args the command line arguments
+   */
+  @SuppressWarnings({ "checkstyle:magicnumber" })
+  public static void main(final String[] args) {
+    // Create a problem instance.
+    int nT = 20;   // # of transmitters
+    int nU = 60;   // # of users
+    int maxC = 7;  // maximum number of clusters
+    int maxS = 5;  // maximum cluster size
+    Problem prob = new Problem(nT, nU, maxC, maxS, 456);
+
+    // Create and solve a MIP model.
+    double timeLimit = 300;  // time limit in seconds for CPLEX
+    try {
+      MIPModel mip = new MIPModel(prob);
+      double obj = mip.solve(timeLimit);
+      System.out.println("\nFinal objective value = " + obj);
+      // Get and summarize the optimal transmitter clusters.
+      Collection<Collection<Integer>> clusters = mip.getClusters();
+      System.out.println("\n" + prob.report(clusters));
+    } catch (IloException ex) {
+      System.out.println("\nCPLEX blew up:\n" + ex.getMessage());
+    }
+
+    // Try the greedy heuristic.
+    Greedy greedy = new Greedy(prob);
+    Collection<Collection<Integer>> g = greedy.getSolution();
+    if (g.isEmpty()) {
+      System.out.println("\nThe greedy heuristic failed.");
+    } else {
+      System.out.println("\nGreedy heuristic results:\n" + prob.report(g));
+    }
+    System.out.println("");
+
+    // Try the revised MIP model.
+//    try {
+//      MIPModel2 mip = new MIPModel2(prob);
+//      double obj = mip.solve(timeLimit);
+//      System.out.println("\nFinal objective value = " + obj);
+//      // Get and summarize the optimal transmitter clusters.
+//      Collection<Collection<Integer>> clusters = mip.getClusters();
+//      System.out.println("\n" + prob.report(clusters));
+//    } catch (IloException ex) {
+//      System.out.println("\nCPLEX blew up:\n" + ex.getMessage());
+//    }
+  }
+
+}