sympy (python) 找到两个方程之间最近的点
sympy (python) Find closest points between 2 equations
我想在 2 个方程上找到 2 个最近的点
我将使用:
python 和同情
我的问题示例:
假设我有 2 个方程式:
y=x^2
y=-(x-3)^2
我想找到这两个方程之间最接近的 2 个点(每个方程一个点)。
我将如何在 sympy 中做到这一点?
注意:我希望能够解决任何n次多项式的问题,例如y=x^4 + x^2 + x
和y=x^3 + 3x^2
等
如果我得到多个解决方案就可以了。
我四处搜索并查看了 https://math.stackexchange.com/questions/731207/closest-distance-between-two-quadratic-curves#:~:text=First%20simply%20take%20the%20derivatives,line%20between%20the%20two%20curves。和其他一些,但我不明白如何以这种方式实现这些解决方案,以便它适用于任何一对 n 次多项式。
我的理解是我必须做一些与垂直于函数切线的线相关的计算,但我不知道如何在 sympy 中做/实现它。
如果这在 sympy 中是不可能的,请指点我一些其他的解决方案,让我找到两个函数之间最近的点。
TLDR:在 sympy 中给定 2 个方程,在 2 个方程上找到 2 个点,其中 2 个点最接近。
谢谢
您可以这样设置要求解的方程式:
In [94]: x1 = symbols('x1')
In [95]: x2 = symbols('x2')
In [96]: y1 = x1**2
In [97]: y2 = -(x2 - 3)**2
In [98]: distance_squared = (y1 - y2)**2 + (x1 - x2)**2
In [99]: distance_squared
Out[99]:
2
2 ⎛ 2 2⎞
(x₁ - x₂) + ⎝x₁ + (x₂ - 3) ⎠
In [100]: eq1, eq2 = [distance_squared.diff(x) for x in [x1, x2]]
In [101]: eq1
Out[101]:
⎛ 2 2⎞
4⋅x₁⋅⎝x₁ + (x₂ - 3) ⎠ + 2⋅x₁ - 2⋅x₂
In [102]: eq2
Out[102]:
⎛ 2 2⎞
-2⋅x₁ + 2⋅x₂ + ⎝x₁ + (x₂ - 3) ⎠⋅(4⋅x₂ - 12)
现在您在 x1
和 x2
中有了一对多项式方程。你只需要找到这些方程的解。请注意,由于 Abel-Ruffini 定理,根据多项式的次数,您可以在此处执行的操作受到限制。可以给出数值解的一般方法是使用 Groebner 基:
In [106]: g = groebner([eq1, eq2], [x1, x2])
In [107]: g
Out[107]:
⎛⎡ 4 3 2 5 4 3
GroebnerBasis⎝⎣181⋅x₁ - 24⋅x₂ + 212⋅x₂ - 624⋅x₂ + 737⋅x₂ - 432, 8⋅x₂ - 96⋅x₂ + 472⋅x₂ - 1206⋅
2 ⎤ ⎞
x₂ + 1656⋅x₂ - 999⎦, x₁, x₂, domain=ℤ, order=lex⎠
In [108]: r = real_roots(g[-1])
In [109]: r
Out[109]:
⎡ ⎛ 3 2 ⎞⎤
⎣CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠⎦
In [112]: solve([g[0], x2 - r[0].n()], [x1, x2]) # numerical result
Out[112]: [(0.728082123067953, 2.27191787693205)]
In [111]: solve([g[0], x2 - r[0]], [x1, x2]) # exact result in terms of rootof
Out[111]:
⎡⎛ 3
⎢⎜ ⎛ 3 2 ⎞ ⎛ 3 2 ⎞
⎢⎜- 212⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ - 737⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ + 432
⎢⎜─────────────────────────────────────────────────────────────────────────────────────────────────
⎣⎝ 181
4 2
⎛ 3 2 ⎞ ⎛ 3 2 ⎞
+ 24⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ + 624⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠
──────────────────────────────────────────────────────────────────────────────────────────, CRootOf
⎞⎤
⎟⎥
⎛ 3 2 ⎞⎟⎥
⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠⎟⎥
⎠⎦
在像这样的简单情况下,可以使用 solve 符号式求解方程,但答案可能非常复杂:
In [114]: solve([eq1, eq2], [x1, x2])
Out[114]:
⎡⎛ __
⎢⎜ ╱ 8
⎢⎜ 3 _____________ 5/6 3 ___________ 3 ╱ ─
⎢⎜- 15066259⋅╲╱ 81 + 3⋅√753 - 549045⋅√251⋅3 ⋅╲╱ 27 + √753 + 410109939 + 14945237⋅√753 ╲╱ 8
⎢⎜───────────────────────────────────────────────────────────────────────────────────────, - ──────
⎢⎜ 2/3 ⎛ 6 ___ 2/3⎞
⎢⎜ (27 + √753) ⋅⎝1647135⋅√251⋅╲╱ 3 + 15066259⋅3 ⎠
⎢⎜
⎣⎝
___________ ⎞ ⎛
1 3⋅√753 ⎟ ⎜
─ + ────── ⎟ ⎜ 3 _____________ 5/6 3 ___________
8 1 ⎟ ⎜1647135⋅√251⋅╲╱ 81 + 3⋅√753 + 15066259⋅3 ⋅╲╱ 27 + √753
─────────── + ─────────────────── + 3⎟, ⎜──────────────────────────────────────────────────────────
3 _____________ ⎟ ⎜ 2/3 ⎛
╱ 81 3⋅√753 ⎟ ⎜ (27 + √753) ⋅⎝1647135⋅√25
2⋅3 ╱ ── + ────── ⎟ ⎜
╲╱ 8 8 ⎠ ⎝
5/6 3 ___________ 3 _____________
+ 549045⋅√251⋅3 ⋅ⅈ⋅╲╱ 27 + √753 + 15066259⋅ⅈ⋅╲╱ 81 + 3⋅√753 + 820219878⋅ⅈ + 29890474⋅√753⋅ⅈ
────────────────────────────────────────────────────────────────────────────────────────────────, 3
2/3 6 ___ 2/3 6 ___ ⎞
1⋅3 + 45198777⋅╲╱ 3 - 15066259⋅3 ⋅ⅈ - 1647135⋅√251⋅╲╱ 3 ⋅ⅈ⎠
_____________ ⎞ ⎛
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753 ⎟ ⎜
⎜- ─ + ────⎟⋅3 ╱ ── + ────── ⎟ ⎜ 3 ____________
⎝ 2 2 ⎠ ╲╱ 8 8 1 ⎟ ⎜1647135⋅√251⋅╲╱ 81 + 3⋅√753
- ────────────────────────────── + ────────────────────────────────⎟, ⎜───────────────────────────
3 _____________⎟ ⎜
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753 ⎟ ⎜
2⋅⎜- ─ + ────⎟⋅3 ╱ ── + ────── ⎟ ⎜
⎝ 2 2 ⎠ ╲╱ 8 8 ⎠ ⎝
_ 5/6 3 ___________ 3 _____________
+ 15066259⋅3 ⋅╲╱ 27 + √753 - 29890474⋅√753⋅ⅈ - 820219878⋅ⅈ - 15066259⋅ⅈ⋅╲╱ 81 + 3⋅√753 - 5490
───────────────────────────────────────────────────────────────────────────────────────────────────
2/3 ⎛ 2/3 6 ___ 6 ___ 2/3 ⎞
(27 + √753) ⋅⎝1647135⋅√251⋅3 + 45198777⋅╲╱ 3 + 1647135⋅√251⋅╲╱ 3 ⋅ⅈ + 15066259⋅3 ⋅ⅈ⎠
_____________
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753
5/6 3 ___________ ⎜- ─ - ────⎟⋅3 ╱ ── + ──────
45⋅√251⋅3 ⋅ⅈ⋅╲╱ 27 + √753 1 ⎝ 2 2 ⎠ ╲╱ 8 8
────────────────────────────, 3 + ──────────────────────────────── - ──────────────────────────────
_____________ 3
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753
2⋅⎜- ─ - ────⎟⋅3 ╱ ── + ──────
⎝ 2 2 ⎠ ╲╱ 8 8
⎞ ⎤
⎟ ⎥
⎟ ⎥
⎟ ⎛3 3⋅ⅈ 3 3⋅ⅈ⎞ ⎛3 3⋅ⅈ 3 3⋅ⅈ⎞⎥
⎟, ⎜─ - ───, ─ - ───⎟, ⎜─ + ───, ─ + ───⎟⎥
⎟ ⎝2 2 2 2 ⎠ ⎝2 2 2 2 ⎠⎥
⎟ ⎥
⎟ ⎥
⎠ ⎦
我想在 2 个方程上找到 2 个最近的点
我将使用:
python 和同情
我的问题示例:
假设我有 2 个方程式:
y=x^2
y=-(x-3)^2
我想找到这两个方程之间最接近的 2 个点(每个方程一个点)。 我将如何在 sympy 中做到这一点?
注意:我希望能够解决任何n次多项式的问题,例如y=x^4 + x^2 + x
和y=x^3 + 3x^2
等
如果我得到多个解决方案就可以了。
我四处搜索并查看了 https://math.stackexchange.com/questions/731207/closest-distance-between-two-quadratic-curves#:~:text=First%20simply%20take%20the%20derivatives,line%20between%20the%20two%20curves。和其他一些,但我不明白如何以这种方式实现这些解决方案,以便它适用于任何一对 n 次多项式。
我的理解是我必须做一些与垂直于函数切线的线相关的计算,但我不知道如何在 sympy 中做/实现它。
如果这在 sympy 中是不可能的,请指点我一些其他的解决方案,让我找到两个函数之间最近的点。
TLDR:在 sympy 中给定 2 个方程,在 2 个方程上找到 2 个点,其中 2 个点最接近。
谢谢
您可以这样设置要求解的方程式:
In [94]: x1 = symbols('x1')
In [95]: x2 = symbols('x2')
In [96]: y1 = x1**2
In [97]: y2 = -(x2 - 3)**2
In [98]: distance_squared = (y1 - y2)**2 + (x1 - x2)**2
In [99]: distance_squared
Out[99]:
2
2 ⎛ 2 2⎞
(x₁ - x₂) + ⎝x₁ + (x₂ - 3) ⎠
In [100]: eq1, eq2 = [distance_squared.diff(x) for x in [x1, x2]]
In [101]: eq1
Out[101]:
⎛ 2 2⎞
4⋅x₁⋅⎝x₁ + (x₂ - 3) ⎠ + 2⋅x₁ - 2⋅x₂
In [102]: eq2
Out[102]:
⎛ 2 2⎞
-2⋅x₁ + 2⋅x₂ + ⎝x₁ + (x₂ - 3) ⎠⋅(4⋅x₂ - 12)
现在您在 x1
和 x2
中有了一对多项式方程。你只需要找到这些方程的解。请注意,由于 Abel-Ruffini 定理,根据多项式的次数,您可以在此处执行的操作受到限制。可以给出数值解的一般方法是使用 Groebner 基:
In [106]: g = groebner([eq1, eq2], [x1, x2])
In [107]: g
Out[107]:
⎛⎡ 4 3 2 5 4 3
GroebnerBasis⎝⎣181⋅x₁ - 24⋅x₂ + 212⋅x₂ - 624⋅x₂ + 737⋅x₂ - 432, 8⋅x₂ - 96⋅x₂ + 472⋅x₂ - 1206⋅
2 ⎤ ⎞
x₂ + 1656⋅x₂ - 999⎦, x₁, x₂, domain=ℤ, order=lex⎠
In [108]: r = real_roots(g[-1])
In [109]: r
Out[109]:
⎡ ⎛ 3 2 ⎞⎤
⎣CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠⎦
In [112]: solve([g[0], x2 - r[0].n()], [x1, x2]) # numerical result
Out[112]: [(0.728082123067953, 2.27191787693205)]
In [111]: solve([g[0], x2 - r[0]], [x1, x2]) # exact result in terms of rootof
Out[111]:
⎡⎛ 3
⎢⎜ ⎛ 3 2 ⎞ ⎛ 3 2 ⎞
⎢⎜- 212⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ - 737⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ + 432
⎢⎜─────────────────────────────────────────────────────────────────────────────────────────────────
⎣⎝ 181
4 2
⎛ 3 2 ⎞ ⎛ 3 2 ⎞
+ 24⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠ + 624⋅CRootOf⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠
──────────────────────────────────────────────────────────────────────────────────────────, CRootOf
⎞⎤
⎟⎥
⎛ 3 2 ⎞⎟⎥
⎝4⋅x - 36⋅x + 110⋅x - 111, 0⎠⎟⎥
⎠⎦
在像这样的简单情况下,可以使用 solve 符号式求解方程,但答案可能非常复杂:
In [114]: solve([eq1, eq2], [x1, x2])
Out[114]:
⎡⎛ __
⎢⎜ ╱ 8
⎢⎜ 3 _____________ 5/6 3 ___________ 3 ╱ ─
⎢⎜- 15066259⋅╲╱ 81 + 3⋅√753 - 549045⋅√251⋅3 ⋅╲╱ 27 + √753 + 410109939 + 14945237⋅√753 ╲╱ 8
⎢⎜───────────────────────────────────────────────────────────────────────────────────────, - ──────
⎢⎜ 2/3 ⎛ 6 ___ 2/3⎞
⎢⎜ (27 + √753) ⋅⎝1647135⋅√251⋅╲╱ 3 + 15066259⋅3 ⎠
⎢⎜
⎣⎝
___________ ⎞ ⎛
1 3⋅√753 ⎟ ⎜
─ + ────── ⎟ ⎜ 3 _____________ 5/6 3 ___________
8 1 ⎟ ⎜1647135⋅√251⋅╲╱ 81 + 3⋅√753 + 15066259⋅3 ⋅╲╱ 27 + √753
─────────── + ─────────────────── + 3⎟, ⎜──────────────────────────────────────────────────────────
3 _____________ ⎟ ⎜ 2/3 ⎛
╱ 81 3⋅√753 ⎟ ⎜ (27 + √753) ⋅⎝1647135⋅√25
2⋅3 ╱ ── + ────── ⎟ ⎜
╲╱ 8 8 ⎠ ⎝
5/6 3 ___________ 3 _____________
+ 549045⋅√251⋅3 ⋅ⅈ⋅╲╱ 27 + √753 + 15066259⋅ⅈ⋅╲╱ 81 + 3⋅√753 + 820219878⋅ⅈ + 29890474⋅√753⋅ⅈ
────────────────────────────────────────────────────────────────────────────────────────────────, 3
2/3 6 ___ 2/3 6 ___ ⎞
1⋅3 + 45198777⋅╲╱ 3 - 15066259⋅3 ⋅ⅈ - 1647135⋅√251⋅╲╱ 3 ⋅ⅈ⎠
_____________ ⎞ ⎛
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753 ⎟ ⎜
⎜- ─ + ────⎟⋅3 ╱ ── + ────── ⎟ ⎜ 3 ____________
⎝ 2 2 ⎠ ╲╱ 8 8 1 ⎟ ⎜1647135⋅√251⋅╲╱ 81 + 3⋅√753
- ────────────────────────────── + ────────────────────────────────⎟, ⎜───────────────────────────
3 _____________⎟ ⎜
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753 ⎟ ⎜
2⋅⎜- ─ + ────⎟⋅3 ╱ ── + ────── ⎟ ⎜
⎝ 2 2 ⎠ ╲╱ 8 8 ⎠ ⎝
_ 5/6 3 ___________ 3 _____________
+ 15066259⋅3 ⋅╲╱ 27 + √753 - 29890474⋅√753⋅ⅈ - 820219878⋅ⅈ - 15066259⋅ⅈ⋅╲╱ 81 + 3⋅√753 - 5490
───────────────────────────────────────────────────────────────────────────────────────────────────
2/3 ⎛ 2/3 6 ___ 6 ___ 2/3 ⎞
(27 + √753) ⋅⎝1647135⋅√251⋅3 + 45198777⋅╲╱ 3 + 1647135⋅√251⋅╲╱ 3 ⋅ⅈ + 15066259⋅3 ⋅ⅈ⎠
_____________
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753
5/6 3 ___________ ⎜- ─ - ────⎟⋅3 ╱ ── + ──────
45⋅√251⋅3 ⋅ⅈ⋅╲╱ 27 + √753 1 ⎝ 2 2 ⎠ ╲╱ 8 8
────────────────────────────, 3 + ──────────────────────────────── - ──────────────────────────────
_____________ 3
⎛ 1 √3⋅ⅈ⎞ ╱ 81 3⋅√753
2⋅⎜- ─ - ────⎟⋅3 ╱ ── + ──────
⎝ 2 2 ⎠ ╲╱ 8 8
⎞ ⎤
⎟ ⎥
⎟ ⎥
⎟ ⎛3 3⋅ⅈ 3 3⋅ⅈ⎞ ⎛3 3⋅ⅈ 3 3⋅ⅈ⎞⎥
⎟, ⎜─ - ───, ─ - ───⎟, ⎜─ + ───, ─ + ───⎟⎥
⎟ ⎝2 2 2 2 ⎠ ⎝2 2 2 2 ⎠⎥
⎟ ⎥
⎟ ⎥
⎠ ⎦