| |||||||||||
| |||||||||||
Python notebooks Description Python notebooks available in the GitHub repository python-notebooks.
Source Files By clicking on a file name, a preview is opened at the bottom of this page. pairwise_distance.ipynb
{
"cells": [
{
"cell_type": "markdown",
"id": "90700428",
"metadata": {},
"source": [
"***pairwise_distance.ipynb***\n",
"\n",
"Determine the positions of $N$ points in $D$-dimensional space as to minimize the ratio of pairwise distances between points, that is, minimize the ratio of *min* and *max* squared distances.\n",
"\n",
"© Copyright 2025 Fair Isaac Corporation\n",
"\n",
"Licensed under the Apache License, Version 2.0 (the \"License\"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0.\n",
" \n",
"Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an \"AS IS\" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.\n",
"\n",
"This example uses FICO® Xpress software. By running it, you agree to the Community License terms of the [Xpress Shrinkwrap License Agreement](https://community.fico.com/s/contentdocument/06980000002h0i5AAA) with respect to the FICO® Xpress software. See the [licensing options](https://www.fico.com/en/fico-xpress-trial-and-licensing-options) overview for additional details and information about obtaining a paid license."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "52204f86",
"metadata": {},
"outputs": [],
"source": [
"# Install the xpress package\n",
"%pip install -q xpress matplotlib plotly"
]
},
{
"cell_type": "markdown",
"id": "68fa499f",
"metadata": {},
"source": [
"## Problem description and formulation"
]
},
{
"cell_type": "markdown",
"id": "5a03e29a",
"metadata": {},
"source": [
"This optimization problem determines the positions of **$N$** points in a **$D$**-dimensional space such that the ratio of the maximum to minimum squared pairwise distances is minimized, which is addressed by **maximizing the minimum distance** in our formulation. This promotes a uniform distribution of points, minimizing clustering and maximizing spatial efficiency.\n",
"\n",
"The **decision variables** are defined below, with $\\mathcal{N}$ and $\\mathcal{D}$ being the corresponding sets of points and dimensions, respectively:\n",
"- $x_{i,k}$: Coordinate of point $i \\in \\mathcal{N}$ in dimension $k \\in \\mathcal{D}$\n",
"- $t_{\\min}$: Minimum squared distance between any two points\n",
"- $t_{\\max}$: Maximum squared distance between any two points (fixed to 1)\n",
"- $r$: Ratio variable, defined as $r = \\frac{1}{t_{\\min}}$\n",
"\n",
"For all pairs of points $i,j \\in \\mathcal{N}$ where $i < j$, the squared Euclidean distance is **constrained** as:\n",
"\n",
"$$ \\sum_{k \\in \\mathcal{D}} (x_{i,k} - x_{j,k})^2 \\geq t_{\\min}, \\qquad \\forall i, j \\in \\mathcal{N}, i < j $$\n",
"$$ \\sum_{k \\in \\mathcal{D}} (x_{i,k} - x_{j,k})^2 \\leq 1, \\qquad \\forall i, j \\in \\mathcal{N}, i < j $$\n",
"\n",
"Additionally, we define the ratio variable by adding the constraint $r \\cdot t_{\\min} = 1$. Since the problem is *scale-invariant*, that is, the outcome of the problem does not change when all lengths are multiplied by the same constant factor, for simplicity we fix $t_{\\max} = 1$.\n",
"\n",
"The **objective function** is to maximize the minimum squared distance: $\\max t_{\\min}$. This is equivalent to minimizing the ratio $\\frac{t_{\\max}}{t_{\\min}}$, since $t_{\\max}$ is fixed to 1. This reformulates the problem as a quadratic problem that can be handled by a specific solver (Xpress-SLP) instead of the generic Global solver. The optimization is subject to a time limit defined by `LIMIT`.\n",
"\n",
"The result is a configuration of $N$ points in $D$-dimensional space that are as evenly spaced as possible, minimizing the ratio of largest to smallest pairwise distances.\n"
]
},
{
"cell_type": "markdown",
"id": "f7e9ef6f",
"metadata": {},
"source": [
"## Model implementation"
]
},
{
"cell_type": "markdown",
"id": "15e0ad67",
"metadata": {},
"source": [
"The previous model is implemented below for the 3D case, with a time limit of 10 seconds for the solver. You can also work on the 2D space, and experiment with different numbers of points and increase or decrease the solver [time limit](https://www.fico.com/fico-xpress-optimization/docs/latest/solver/optimizer/HTML/TIMELIMIT.html). A lower bound of 0.1 is defined for $t_{min}$ to ensure that points are sufficiently far from each other, giving a maximum ratio of approximately $10$ for $N = 12$ points.\n",
"\n",
"After the optimization, the objective (ratio) and solution values are printed following the solver logs (ON by default)."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "ff48ce99",
"metadata": {},
"outputs": [],
"source": [
"import xpress as xp\n",
"import matplotlib.pyplot as plt\n",
"import plotly.graph_objects as go\n",
"import plotly.io as pio\n",
"\n",
"N = 12 # Number of points to be placed.\n",
"D = 3 # Number of dimensions: 3 for 3D, 2 for 2D.\n",
"LIMIT = 10 # Time limit for the solver.\n",
"\n",
"POINTS = range(N) # Points to be placed.\n",
"DIMS = range(D) # Dimensions.\n",
"\n",
"# Create a problem instance.\n",
"p = xp.problem()\n",
"\n",
"# Decision variables.\n",
"x = p.addVariables(POINTS, DIMS, ub=1, name=\"x\") # N points in d-dimensional space\n",
"t_min = p.addVariable(lb=0.1, name=\"t_min\") # Minimum squared distance variable.\n",
"t_max = p.addVariable(name=\"t_max\") # Maximum squared distance variable.\n",
"r = p.addVariable(name=\"r\") # Ratio variable.\n",
"\n",
"# Constraints to ensure that the points are placed in a way that\n",
"# the minimum distance between any two points is t_min and the maximum is 1\n",
"p.addConstraint(xp.Sum((x[i,k] - x[j,k])**2 for k in DIMS) >= t_min for i in POINTS for j in POINTS if i < j)\n",
"p.addConstraint(xp.Sum((x[i,k] - x[j,k])**2 for k in DIMS) <= 1 for i in POINTS for j in POINTS if i < j)\n",
"\n",
"# We wish to minimize t_max/t_min, which is invariant for scaling. We can\n",
"# thus narrow the search to configurations where t_max is 1\n",
"p.addConstraint(t_max == 1)\n",
"# Compute ratio\n",
"p.addConstraint(r * t_min == 1)\n",
"\n",
"# Objective: minimize the square root of ratio t_max / t_min.\n",
"# Instead of minimizing r=1/t_min, we can improve the formulation by\n",
"# maximizing t_min.\n",
"p.setObjective(t_min, sense=xp.maximize)\n",
"\n",
"# Set the solver parameters.\n",
"p.controls.timelimit = LIMIT\n",
"\n",
"# Solve the problem.\n",
"p.optimize()\n",
"\n",
"# Print the solution\n",
"status = p.attributes.solstatus\n",
"if status in {xp.SolStatus.FEASIBLE, xp.SolStatus.OPTIMAL}:\n",
" print(\"Solution value (ratio):\", p.getSolution(r))\n",
" print(\"t_max:\", p.getSolution(t_max), \", t_min:\", p.getSolution(t_min))\n",
" xsol = p.getSolution(x)\n",
" for i in POINTS:\n",
" for j in DIMS:\n",
" print(f\"{x[i,j].name}: {xsol[i,j]}\")\n",
"else:\n",
" print(f\"No solution found with time limit of {LIMIT} sec\")\n"
]
},
{
"cell_type": "markdown",
"id": "bfa4ee35",
"metadata": {},
"source": [
"## Visualization"
]
},
{
"cell_type": "markdown",
"id": "6ceb030f",
"metadata": {},
"source": [
"The code cell below generates a plot in either 2D or 3D (depending on the value of $D$ chosen above) using *matplotlib*. \n",
"\n",
"The plot containts the solution points, and the lines between any two points placed at the **minimum distance (blue)** or the **maximum distance (red)** amongst the whole set of pairwise distances."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "90bb8e2b",
"metadata": {},
"outputs": [],
"source": [
"print(f\"N = {N}, D = {D}, Ratio = {p.getSolution(r)}\")\n",
"\n",
"# Convert to dictionary with point ID as key and list of coordinates as value\n",
"points = {}\n",
"for (point_id, dimension), value in xsol.items():\n",
" if point_id not in points:\n",
" # Initialize list with None values up to the maximum dimension\n",
" max_dim = max(d for _, d in xsol.keys()) + 1\n",
" points[point_id] = [None] * max_dim\n",
" points[point_id][dimension] = value\n",
"\n",
"# Calculate squared distance\n",
"def squared_distance(p1, p2, ndims):\n",
" return sum((p1[i] - p2[i])**2 for i in range(ndims))\n",
"\n",
"# Create a plot\n",
"fig = plt.figure(figsize=(10, 6))\n",
"ax = fig.add_subplot(111, projection='3d' if D == 3 else None)\n",
"\n",
"# Plot points\n",
"for point in points.values():\n",
" ax.scatter(*point, color='black')\n",
"\n",
"# Draw lines with color based on distance\n",
"tol = p.controls.feastol\n",
"for i in points:\n",
" for j in points:\n",
" p1, p2 = points[i], points[j]\n",
" dist = squared_distance(p1, p2, D)\n",
" if dist - p.getSolution(t_min) <= tol:\n",
" if D == 2: # For 2D, we can plot in 2D space\n",
" ax.plot([p1[0], p2[0]], [p1[1], p2[1]], color='blue')\n",
" else: # For 3D, we plot in 3D space\n",
" ax.plot([p1[0], p2[0]], [p1[1], p2[1]], [p1[2], p2[2]], color='blue')\n",
" if 1 - dist <= tol: \n",
" if D == 2:\n",
" ax.plot([p1[0], p2[0]], [p1[1], p2[1]], color='red')\n",
" else:\n",
" ax.plot([p1[0], p2[0]], [p1[1], p2[1]], [p1[2], p2[2]], color='red')\n",
"\n",
"ax.set_xlabel('x')\n",
"ax.set_ylabel('y')\n",
"ax.set_zlabel('z') if D == 3 else None\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "cb6b8f1d",
"metadata": {},
"source": [
"## Interactive 3D visualization"
]
},
{
"cell_type": "markdown",
"id": "1e20b8ea",
"metadata": {},
"source": [
"If you selected $D = 3$, then you can generate an **interactive 3D plot** of the same solution with the code cell below, which uses the *plotly* library. By default, the renderer is setup for Visual Studio Code, but you can adapt this to your own environment.\n",
"\n",
"After rendering, you can do *orbital rotation of the plot in any direction* and see the solution from different angles, as well as observing the coordinates of each point when hovering over. An `HMTL` file named `interactive_plot.html` is also exported in case you prefer to explore it in a web browser.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "79556195",
"metadata": {},
"outputs": [],
"source": [
"# Set the default renderer for Plotly to 'vscode'\n",
"pio.renderers.default = 'vscode' # 'notebook', 'iframe', 'browser', ...\n",
"\n",
"# Calculate squared distance\n",
"def squared_distance(p1, p2):\n",
" return sum((p1[i] - p2[i])**2 for i in range(len(p1)))\n",
"\n",
"# Create a 3D scatter plot\n",
"fig = go.Figure()\n",
"\n",
"# Add points to the plot\n",
"for point in points.values():\n",
" fig.add_trace(go.Scatter3d(\n",
" x=[point[0]], y=[point[1]], z=[point[2]],\n",
" mode='markers',\n",
" marker=dict(color='black', size=5)\n",
" ))\n",
"\n",
"# Draw lines with color based on distance\n",
"tol = p.controls.feastol\n",
"for i in points:\n",
" for j in points:\n",
" if i != j:\n",
" p1, p2 = points[i], points[j]\n",
" dist = squared_distance(p1, p2)\n",
" if dist - p.getSolution(t_min) <= tol:\n",
" fig.add_trace(go.Scatter3d(\n",
" x=[p1[0], p2[0]], y=[p1[1], p2[1]], z=[p1[2], p2[2]],\n",
" mode='lines',\n",
" line=dict(color='blue')\n",
" ))\n",
" if 1 - dist <= tol:\n",
" fig.add_trace(go.Scatter3d(\n",
" x=[p1[0], p2[0]], y=[p1[1], p2[1]], z=[p1[2], p2[2]],\n",
" mode='lines',\n",
" line=dict(color='red')\n",
" ))\n",
"\n",
"# Set axis labels\n",
"fig.update_layout(scene=dict(\n",
" xaxis_title='x',\n",
" yaxis_title='y',\n",
" zaxis_title='z'),\n",
" showlegend=False, # Remove legend entries\n",
" width=1000, # Set figure width\n",
" height=800 # Set figure height\n",
")\n",
"\n",
"# Save the plot as an HTML file\n",
"fig.write_html(\"interactive_plot.html\")\n",
"\n",
"# Show the plot\n",
"fig.show()\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.12.9"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
| |||||||||||
| © Copyright 2025 Fair Isaac Corporation. |
