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()

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.