求解参数曲线之间的交点

Question

两条曲线的参数方程如下：

Curve1: r(t) = (2(t-sin(t)),2(1 -cos(t)))

Curve2: s(t) = (2t - sin(t),2 - cos(t))

我需要找到区域[0,4π]中的相交点。

我能够绘制上述区域的图形，并观察到4个交点。但我无法确定确切的交点。

对于非参数方程，可以使用sympy中的fsolve，但以其参数形式给出的曲线，我无法找到解决方法。

t = np.arange(-0.25*np.pi,4.25*np.pi,0.01)
rx = np.zeros(len(t))
ry = np.zeros(len(t))
for i in t:
rx = 2*(t - np.sin(t))
ry = 2*(1 - np.cos(t))
sx = np.zeros(len(t))
sy = np.zeros(len(t))
for i in t:
sx = 2*t - np.sin(t)
sy = 2 - np.cos(t)
plt.plot(rx,ry)
plt.plot(sx,sy)

Answer 1

对于给定的x，您可以为每条曲线找到t，并查看相应的y是否相同。你可以用一些网格在x范围内寻找e曲线碰撞的位置，并使用平分法在更精确的x上归零。由于无法求解t的参数x(t) - x，因此必须使用nsolve来找到近似的t。在将Curve1的OP方程校正为与您随后编写的代码中的值相同后，类似这样的操作会找到您的4个根的值(以图形方式确认(。

f = lambda xx: a[1].subs(t, tt)-b[1].subs(t,nsolve(b[0]-xx,tlast))
tlast = 0 # guess for t for a given xx to be updated as we go
tol = 1e-9  # how tight the bounds on x must be for a solution
dx = 0.1
for ix in range(300):
xx = ix*dx
tt=nsolve(a[0]-xx,tlast)
y2 = f(xx)
if ix != 0 and yold*y2 < 0 and tt<4*pi:
tlast = tt  # updating guess for t
# bisect for better xx now that bounding xx are found
x1 = xx-dx
x2 = xx
y1 = yold
while x2 - x1 > tol:
xm = (x1 + x2)/2
ym = f(xm)
if ym*y1 < 0:
y2 = ym
x2 = xm
elif ym != 0:
y1 = ym
x1 = xm
else:
break
print(xm)  # a solution
yold = y2

我不知道有什么比SymPy更自动化的方法可以做到这一点。

Answer 2

两条曲线不同时相交(这将是sin(t(=cos(t(=0的点，没有解(。所以你真的想知道什么时候

R = (2*t1 - 2*sin(t1), 2 - 2*cos(t1))
S = (2*t2 - sin(t2), 2 - cos(t2))

相交。

这是两个有两个未知数的方程，所以用sympy.nsolve求解很简单。你必须稍微摆弄一下起始值，才能找到收敛到不同解决方案的值。如果你从图表中大致知道它们是什么，那就是最好的开始。

>>> t1, t2 = symbols('t1 t2')
>>> R = (2*t1 - 2*sin(t1), 2 - 2*cos(t1))
>>> S = (2*t2 - sin(t2), 2 - cos(t2))
>>> nsolve([R[0] - S[0], R[1] - S[1]], [t1, t2], [1, 1])
Matrix([
[ 1.09182358380672],
[0.398264297579454]])
>>> nsolve([R[0] - S[0], R[1] - S[1]], [t1, t2], [5, 5])
Matrix([
[5.19136172337286],
[5.88492100960013]])
>>> nsolve([R[0] - S[0], R[1] - S[1]], [t1, t2], [7, 7])
Matrix([
[7.37500889098631],
[6.68144960475904]])
>>> nsolve([R[0] - S[0], R[1] - S[1]], [t1, t2], [10, 10])
Matrix([
[11.4745470305524],
[12.1681063167797]])