【Python】ランダム配置の円に対する接点の抽出

Question

前提・実現したいこと

現在、Pythonを用いた円のRandom Packingに取り組んでいます。
現時点では、指定個数かつ半径が異なる円を重ならないように配置するプログラムまでは書くことはできました。
下記にコードを示します。

import cvxpy as cvx
import numpy as np
import matplotlib.pyplot as plt
import dccp

n = 40
r = np.linspace(1,5,n)

c = cvx.Variable((n,2))
constr = []
for i in range(n-1):
    for j in range(i+1,n):
        constr.append(cvx.norm(c[i,:]-c[j,:])>=r[i]+r[j])
prob = cvx.Problem(cvx.Minimize(cvx.max(cvx.norm(c, axis=1)+r)), constr)
# prob = cvx.Problem(cvx.Minimize(cvx.max(cvx.max(cvx.abs(c),axis=1)+r)), constr)
prob.solve(method = 'dccp', ccp_times = 1)

l = cvx.max(cvx.max(cvx.abs(c),axis=1)+r).value*2
pi = np.pi
ratio = pi*cvx.sum(cvx.square(r)).value/cvx.square(l).value
print("ratio =", ratio)
print(prob.status)

次の画像はmatplotlibで結果を描画したものです。

この図から以下の情報を抽出したいと考えています。

• 接点の位置
• 2円の接点と中心を通る線の方向ベクトル

なお上記のコードで求めた円は厳密には接していないため、ここでの「接点」は(2円の中心間距離)-(2円の半径の和)<(ある閾値)を満たすものとしたいです。

発生している問題

現在考えた手法としては以下のようなものがあります

1. 各円とそれ以外の円の中心を結ぶ直線の方程式を抽出
2. 1で抽出した直線の方程式に対して、2円を除くすべての円について交差判定を実施
3. 交差しなかった直線のみから方向ベクトルと接点位置を計算

しかしながらこの方法の問題点は計算量がO(n^3)と大きくなってしまうことです。

質問

以上を踏まえまして、下記が質問となります。

1. 上記の計算量を減らすアルゴリズムなどは存在しますでしょうか？
2. １．につきまして参考文献などがありましたら、教えていただけますでしょうか？

追記

今回は円の直径が異なるため、閾値よりも小さい半径を持つ円が存在した場合、(2円の中心間距離)-(2円の半径の和)<(ある閾値)のみで判定を行うと、図のような間違った結果が出てきてしまいます。

これを避けるために交差判定が必要と考え、上記の手法に加えております。。

何卒よろしくお願いいたします。

Answer 1

２個の円が

(2円の中心間距離)-(2円の半径の和)<(ある閾値)

を満たすなら「接点」を計算する……という話ではダメなのでしょうか？
（総当りでも O(n^2) に思える）

Comments

[Kashi_wamochi]

ご回答ありがとうございます。
ご指摘の方法で試してみたのですが、半径が等しい円の場合はその方法で問題ありませんでした。
ただ今回は半径にばらつきがあり、図のような配置の場合、円3も間違って接点として検出されてしまいます...(図は質問内に追記させていただきました)
そのため交差判定が必要なのではと思い、上記のような考えに至りました。

[fana]

> 質問内に追記
I can't find it.

[Kashi_wamochi]

すみません。ただいま追記させていただきました。

[Kashi_wamochi]

追記した図では円1と円0がかなり離れているように見えるかもしれませんが、これは図化のために若干誇張して書かせていただきました。
実際の解析では円1がより大きく、円2や円3がより小さいものになる予定です。

[fana]

正直いまいち話が謎なのですが，
この世界における「"接している"とは何なのか？」というところが不明瞭であるということに思えます．

円のサイズにものすごい開きがある，という話なのであれば，
閾値を定数ではなく，円のサイズに応じて（判定対象たる２つの円の半径に基づいて）可変させるような話で多少は逃げられるようにも思えますが…　果たして？

[fana]

あと，「交差判定」の話とのつながりもわからないです．
（私が，そこの話の意味を理解できていないのだと思いますが…）
例えばその図で，円0と3との関係は
・円2が存在すれば，接触していない
・円2が存在しないならば，接触している
として扱うのでしょうか？（0と3の距離は同じなのに，扱いが変わる？）

[fana]

そういう話（複数の細かい円の並びが，「閾値」が作る帯の幅の中に入り込んでしまうことに問題がある）である場合，
最初にボロノイ図みたいなのを作って，判定対象自体を絞るような話にすれば良いのかもしれません．

（あるいは，許されるならば，円の配置自体をボロノイ図から出発するとか？）

[Kashi_wamochi]

詳しい回答をありがとうございます。
>この世界における「"接している"とは何なのか？」というところが不明瞭であるということに思えます．

こちらの書き方が不適切でした。より正確な書き方をすれば、今回行いたいのは、「接するあるいは隣り合う２円の中心点を結んだ位置ベクトルを求めたい」ということです。ここでの「接するあるいは隣り合う」というのは、下記の条件を満たすものになります。

1. (2円の中心間距離)-(2円の半径の和)<(ある閾値)
2. 中心同士を結んだ線が他の円と交差しない

>あと，「交差判定」の話とのつながりもわからないです．
上の条件のうち、2が「交差判定」と呼ばれるものになります。(勝手に名前を付けておいて、説明をつけておりませんでした。申し訳ありません。)

>最初にボロノイ図みたいなのを作って，判定対象自体を絞るような話にすれば良いのかもしれません．
ありがとうございます。恥ずかしながらボロノイ図について発想が回りませんでした。
試してみたいと思います。

[KAZENOMACHI]

すみませんが、作画をmatplotlibで行ったとありますが、どのようなコードになるのでしょうか？