使用迭代的数值方法查找圆的圆切线
Find a circle circle tangent using an iterative, numerical, method
给定:两个圆的圆心和半径
我的问题是我不知道如何使用数值方法求圆的圆切线。使用 http://mathworld.wolfram.com/Circle-CircleTangents.html 我会知道如何手动解决它,但不知道如何使用 matlab/octave 来解决。
是设置问题让我很为难。即创建方程组。在上面的网站上,他们使用了这样一个事实,即从中心到切点的线之间的点积垂直于切线。但是,我不知道如何将其转化为我的模型。如果我使用matlab的"dot"函数,我将如何设置雅可比矩阵?
使用 MATLAB 内置函数
最简单的方法是按原样定义函数,然后使用 fsolve()
函数。例如(代码未测试):
function y = cctang(t)
%t is a 1-by-4 line vector (first two members = t1, last two members = t2)
%y is a 1-by-4 vector
%Initialization of centers and radii
x1 = [0 0]; r1=4;
x2 = [1 1]; r2=0.5;
y= zeros(1, 4);
y(1) = ( t(3:4) - x2) * ( t(3:4) - t(1:2))';
y(2) = ( t(1:2) - x1) * ( t(3:4) - t(1:2))';
y(3) = norm(t(1:2)-x1)^2 - r1^2;
y(4) = norm(t(3:4)-x2)^2 - r2^2;
end
然后执行:
t0 = [0 0 1 1];
solution = fsolve(@cctang,t0);
使用自定义迭代函数
如果您想从头开始构建所有内容,您可以构建一个自定义雅可比函数并将其输入到您在问题中给出的求解函数中。雅可比函数应如下所示:
function J = myJacobian(t)
%returns a 4-by-4 matrix
x1 = [0 0]; r1=4;
x2 = [1 1]; r2=0.5;
J = zeros(4,4);
J(1,:) = [ x2(1)-t(3), x2(2)-t(4), 2*t(3)-x2(1), 2*t(4)-x2(2)];
%etc...
end
使用符号工具箱
前面方法中的大部分 "hard work" 来自于您需要手动计算导数的事实。
使用 Symbolic Toolbox,您可以用符号定义函数,然后使用 built-in jacobian()
function 计算 Jacobian 矩阵。
这些将是您自定义数值求解函数的输入。在每次迭代中,您应该用实际数值替换符号变量(参见 subs()
)。
这将使您的代码更加通用并且对人为错误不那么敏感(例如,我在上次计算的 one 导数中出错的概率是多少代码块?)。
作为 Yellow 解决方案的替代方案,我想指出该问题可以通过更具几何风格的迭代来解决。此解决方案适用于任何平滑曲线,而不仅仅是圆。
我们称圆心为 A
和 B
,半径为 R[A]
和 R[B]
的圆 C[A]
和 C[B]
。方括号表示下标。
您还有一些线端点的初步估计,可能来自用户选择。调用这些点 P[A][0]
和 P[B][0]
.
对于 P[A][i]
和 P[B][i]
,计算 P[A][i+1]
和 P[B][i+1]
如下:
- 通过将
P[A][i]
和 P[B][i]
投影到圆上来创建直线 L[i]
。
- 让
Q[A][i]
成为 L[i]
上最接近 A
的点。
- 让
Q[B][i]
同样是 L[i]
上最接近 B
的点。
- 设
P[A][i+1]
为Q[A][i]
在C[A]
上的投影。
- 设
P[B][i+1]
为 Q[B][i]
到 C[B]
的投影。
- 如果
||P[A][i+1] - Q[A][i]|| + ||P[A][i+1] - Q[A][i]|| < epsilon
,线段P[A][i+1], P[B][i+1]
就是你的解决方案。
"Deflate" 两个圆直到一个变成一个点。放气时,外切线方向一致
半径为R
的圆,距离圆心D
的一点的两条切线,求解一个直角三角形就够了。特别是,接触点与中心线的位置角由 sin(phi)=R/D
.
给出
通过类似地膨胀较大的圆来找到两条内切线。
更新:简单的解析解。
不失一般性,让缩小的圆以原点为中心,半径为单位(应用平移和缩放)。让 (U, V)
外点和 (X, Y)
圆上的一个点(这样 X²+Y²=1
)。
我们用点积表示切线与半径向量的正交性,U . (U - X) + Y . (V - Y) = 0, i.e. U . X + V . Y = 1
。
让我们通过除以 √(U² + V²)
来归一化得到 u . X + v . Y = t
,我们重写 (v . Y)² = (t - u . X)²
,或 X² - 2 u . t . X + t² - v² = 0
.
解决方法很简单
X = u . t +/- v . √(1 - t²)
Y = v . t -/+ u . √(1 - t²)
乘以两个原始圆的半径,得到从圆心到接触点的向量。
给定:两个圆的圆心和半径
我的问题是我不知道如何使用数值方法求圆的圆切线。使用 http://mathworld.wolfram.com/Circle-CircleTangents.html 我会知道如何手动解决它,但不知道如何使用 matlab/octave 来解决。
是设置问题让我很为难。即创建方程组。在上面的网站上,他们使用了这样一个事实,即从中心到切点的线之间的点积垂直于切线。但是,我不知道如何将其转化为我的模型。如果我使用matlab的"dot"函数,我将如何设置雅可比矩阵?
使用 MATLAB 内置函数
最简单的方法是按原样定义函数,然后使用 fsolve()
函数。例如(代码未测试):
function y = cctang(t)
%t is a 1-by-4 line vector (first two members = t1, last two members = t2)
%y is a 1-by-4 vector
%Initialization of centers and radii
x1 = [0 0]; r1=4;
x2 = [1 1]; r2=0.5;
y= zeros(1, 4);
y(1) = ( t(3:4) - x2) * ( t(3:4) - t(1:2))';
y(2) = ( t(1:2) - x1) * ( t(3:4) - t(1:2))';
y(3) = norm(t(1:2)-x1)^2 - r1^2;
y(4) = norm(t(3:4)-x2)^2 - r2^2;
end
然后执行:
t0 = [0 0 1 1];
solution = fsolve(@cctang,t0);
使用自定义迭代函数
如果您想从头开始构建所有内容,您可以构建一个自定义雅可比函数并将其输入到您在问题中给出的求解函数中。雅可比函数应如下所示:
function J = myJacobian(t)
%returns a 4-by-4 matrix
x1 = [0 0]; r1=4;
x2 = [1 1]; r2=0.5;
J = zeros(4,4);
J(1,:) = [ x2(1)-t(3), x2(2)-t(4), 2*t(3)-x2(1), 2*t(4)-x2(2)];
%etc...
end
使用符号工具箱
前面方法中的大部分 "hard work" 来自于您需要手动计算导数的事实。
使用 Symbolic Toolbox,您可以用符号定义函数,然后使用 built-in jacobian()
function 计算 Jacobian 矩阵。
这些将是您自定义数值求解函数的输入。在每次迭代中,您应该用实际数值替换符号变量(参见 subs()
)。
这将使您的代码更加通用并且对人为错误不那么敏感(例如,我在上次计算的 one 导数中出错的概率是多少代码块?)。
作为 Yellow 解决方案的替代方案,我想指出该问题可以通过更具几何风格的迭代来解决。此解决方案适用于任何平滑曲线,而不仅仅是圆。
我们称圆心为 A
和 B
,半径为 R[A]
和 R[B]
的圆 C[A]
和 C[B]
。方括号表示下标。
您还有一些线端点的初步估计,可能来自用户选择。调用这些点 P[A][0]
和 P[B][0]
.
对于 P[A][i]
和 P[B][i]
,计算 P[A][i+1]
和 P[B][i+1]
如下:
- 通过将
P[A][i]
和P[B][i]
投影到圆上来创建直线L[i]
。 - 让
Q[A][i]
成为L[i]
上最接近A
的点。 - 让
Q[B][i]
同样是L[i]
上最接近B
的点。 - 设
P[A][i+1]
为Q[A][i]
在C[A]
上的投影。 - 设
P[B][i+1]
为Q[B][i]
到C[B]
的投影。 - 如果
||P[A][i+1] - Q[A][i]|| + ||P[A][i+1] - Q[A][i]|| < epsilon
,线段P[A][i+1], P[B][i+1]
就是你的解决方案。
"Deflate" 两个圆直到一个变成一个点。放气时,外切线方向一致
半径为R
的圆,距离圆心D
的一点的两条切线,求解一个直角三角形就够了。特别是,接触点与中心线的位置角由 sin(phi)=R/D
.
通过类似地膨胀较大的圆来找到两条内切线。
更新:简单的解析解。
不失一般性,让缩小的圆以原点为中心,半径为单位(应用平移和缩放)。让 (U, V)
外点和 (X, Y)
圆上的一个点(这样 X²+Y²=1
)。
我们用点积表示切线与半径向量的正交性,U . (U - X) + Y . (V - Y) = 0, i.e. U . X + V . Y = 1
。
让我们通过除以 √(U² + V²)
来归一化得到 u . X + v . Y = t
,我们重写 (v . Y)² = (t - u . X)²
,或 X² - 2 u . t . X + t² - v² = 0
.
解决方法很简单
X = u . t +/- v . √(1 - t²)
Y = v . t -/+ u . √(1 - t²)
乘以两个原始圆的半径,得到从圆心到接触点的向量。