How to create a modified voronoi algorithm for random points with physical restriction

Question

Voronoi algorithm has no doubt provided a amenable approach to divide a plane into regions based on distance to points in a specific subset of the plane. Such the Voronoi diagram of a set of points is dual to its Delaunay triangulation.This goal now can be directly achieved by using the module scipy as

import scipy.spatial

point_coordinate_array = np.array(point_coordinates)
delaunay_mesh = scipy.spatial.Delaunay(point_coordinate_array)
voronoi_diagram = scipy.spatial.Voronoi(point_coordinate_array)

# plot(delaunay_mesh and voronoi_diagram using matplotlib)

when the given points beforehand. The results can be shown in Fig 1,

in which the areas bounded by green dashed lines are the delaunay triangles of all points and the closed region into blue solid lines is of course the voronoi cell of the center point(for a better Visualization, only the closed region is shown here)

Up to now, all things are seen to be perfect. But for actual application, all points may have their own physical meanings. (For example, these points may have 'radii' varible when they represent the natural particles). And the common voronoi algorithm above is to be more or less inappropriate for such the case in which the complex physical restriction may be considered. As shown in Fig 2, the ridges of voronoi cell may intersect the particle's boundary. It can not meet the physical requirment any more.

My question now is how to create a modified voronoi algorithm (maybe it can not called voronoi any more) to deal with this physical restriction. This purpose is shown roughly in Fig 3, the region closed by blue dashed lines is just what I want.

All the pip requirments are:

1.numpy-1.13.3+mkl-cp36-cp36m-win_amd64.whl

2.scipy-0.19.1-cp36-cp36m-win_amd64.whl

3.matplotlib-2.1.0-cp36-cp36m-win_amd64.whl

and all of them can be directly downloaded in http://www.lfd.uci.edu/~gohlke/pythonlibs/

My codes are updated for a better modifcation, they are as

import numpy as np
import scipy.spatial
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
from matplotlib.collections import PatchCollection

# give the point-coordinate array for contribute the tri-network.
point_coordinate_array = np.array([[0,0.5],[8**0.5,8**0.5+0.5],[0,-3.5],[-np.sqrt(15),1.5]])

# give the physical restriction (radius array) here.
point_radius_array = np.array([2.5,1.0,1.0,1.0])

# create the delaunay tri-mesh and voronoi diagram for the given points here.
point_trimesh = scipy.spatial.Delaunay(point_coordinate_array)
point_voronoi = scipy.spatial.Voronoi(point_coordinate_array)

# show the results using matplotlib.
# do the matplotlib setting here.
fig_width = 8.0; fig_length = 8.0
mpl.rc('figure', figsize=((fig_width * 0.3937), (fig_length * 0.3937)), dpi=300)
mpl.rc('axes', linewidth=0.0, edgecolor='red', labelsize=7.5, labelcolor='black', grid=0)
mpl.rc('xtick.major', size=0.0, width=0.0, pad=0)
mpl.rc('xtick.minor', size=0.0, width=0.0, pad=0)
mpl.rc('ytick.major', size=0.0, width=0.0, pad=0)
mpl.rc('ytick.minor', size=0.0, width=0.0, pad=0)
mpl.rc('figure.subplot', left=0.0, right=1.0, bottom=0.065, top=0.995)
mpl.rc('savefig', dpi=300)
ax_1 = plt.figure().add_subplot(1, 1, 1)
plt.gca().set_aspect('equal')
ax_1.set_xlim(-5.5, 8.5)
ax_1.set_ylim(-4.5, 7.5)
ax_1.set_xticklabels([])
ax_1.set_yticklabels([])

# plot all the given points and vertices here.
ax_1.scatter(point_coordinate_array[:,0],point_coordinate_array[:,1],
             s=7.0,c='black')
ax_1.scatter(point_voronoi.vertices[:,0],point_voronoi.vertices[:,1],
             s=7.0,c='blue')

# plot the delaunay tri-mesh here.
ax_1.triplot(point_trimesh.points[:,0],point_trimesh.points[:,1],
             point_trimesh.vertices,
             linestyle='--',dashes=[2.0]*4,color='green',lw=0.5)

# plot the voronoi cell here.(only the closed one)
ax_1.plot(point_voronoi.vertices[:,0],point_voronoi.vertices[:,1],
          lw=1.0,color='blue')
ax_1.plot([point_voronoi.vertices[-1][0],point_voronoi.vertices[0][0]],
          [point_voronoi.vertices[-1][1],point_voronoi.vertices[0][1]],
          lw=1.0,color='blue')

# plot all the particles here.(point+radius)
patches1 = [Circle(point_coordinate_array[i], point_radius_array[i]) 
            for i in range(len(point_radius_array))]
ax_1.add_collection(PatchCollection(patches1, linewidths=1.0, 
                    edgecolor='black', facecolors='none', alpha=1.0))
# save the .png file.
plt.savefig('Fig_a.png',dpi=300)
plt.close()

Bigham

This is a cool question, I have an algorithm written down. Would you post the pip requirements and the plot function so that I can attempt to reproduce your diagrams. – Externality

I think the first step would be to create a Delaunay triangulation, then your first boundary would be a line perpendicular to a line joining the circle centres & going through the point that was mid-way between the circumferences of the 2 circles. Interesting problem. – Feathery

This looks like a weighted Voronoi diagram, or more specifically, a Power Diagram. – Radiosensitive

@Splatmistro Thanks for replying. The pip requirements and the updated codes are provided in edited version of my question. – Bigham

While I do think that this might be an interesting and useful problem, I have problems understanding what the desired output should be. Is Fig3 the desired output? I guess not? But then, what represents the limitation? – Kike

If you define boundaries of Voronoi diagram cells in your case as the set of points for which two closest circles are equidistant, then I believe the diagram is not constructed from lines. You can see it for a point and a circle f.e. So the already mentioned power diagrams (using not the distance to the circle but length of tangent length, and featuring straight cell borders) may be an exclusive option actually. – Professionalism

This might be solved by a Power Diagram. You can look at this question that gives the code to build a Power Diagram in 2D. – Schuh

Answer 1

It should be solved by a Power Diagram approach. You can find algorithms for 2D Power diagrams:

• In a former question
• In GitHub

The Power diagram is constructed from a finite set of circles in the space. The cell for a given circle consists of all the points for which the power distance to the circle is smaller than the power distance to the other circles. The power distance of a point is defined as the square of the point distance to circle center minus the square of circle radius.

Power Diagram have straight cell edges. This seems to correspond to the question. Weighted Voronoi diagrams also exist but can have as edges circle arcs.

Schuh