椭圆形计算的算法

Question

我有n个圆，必须完美地围绕一个椭圆，如下图所示： 在这张图片中，我需要找出每个圆在椭圆周围的位置，并且能够计算出适合这些圆的完美内部椭圆。

我知道的信息是每个圆的半径（都相同）和圆的数量。

希望这次发布的内容更清晰了。谢谢您的帮助。如果需要更多解释，请告诉我。

- user2360915

---

这是一些很棒的ASCII艺术！要发布图片，需将其托管在另一个网站上，并使用Markdown引用它。请参阅此参考资料。 - beatgammit

谢谢。我没有考虑过在其他地方托管图片。我更新了帖子，现在应该更好了。 - user2360915

嗨，Spektre。是的，它必须是一个椭圆。如果它是一个圆，那么我知道使用余弦定律和正多边形的解决方案。如果它是一个椭圆，似乎更困难。必须找到椭圆的参数，并基于“N”个周围圆的半径“r”。 - user2360915

我担心这更像是一个数学问题。仍然对这个问题感兴趣。 - Fabinout

我认为数学是必要的，但由于它没有代数解（我们无法精确计算椭圆的长度或边界多边形），编程是唯一剩下的解决方案（生成/测试，或者在我的情况下尝试添加以提高准确性迭代）。当然，我可能会错过其中的一些基本公式...所以如果有，请与我们分享，让我们也受益。 :) - Spektre

Answer 1

好的，我理解您已经知道圆的公共半径R0和它们的数量N，并想要了解椭圆内部的参数和位置。

如果我们将椭圆转换为圆，则得到以下内容：

const int N=12; // number of satelite circles
const double R=10.0;    // radius of satelite circles
struct _circle { double x,y,r; } circle[N]; // satelite circles

int i;
double x,y,r,l,a,da;
x=0.0;  // start pos of first satelite circle
y=0.0;
r=R;
l=r+r;  // distance ang angle between satelite circle centers
a=0.0*deg;
da=divide(360.0*deg,N);
for (i=0;i<N;i++)
    {
    circle[i].x=x; x+=l*cos(a);
    circle[i].y=y; y+=l*sin(a);
    circle[i].r=r; a+=da;
    }
// inside circle params
_circle c;
r=divide(0.5*l,sin(0.5*da))-R;
c.x=circle[i].x;
c.y=circle[i].y+R+r;
c.r=r;

[编辑1]

对于椭圆，这是一个全新的挑战（花了我两个小时才找到所有的怪癖）

const int    N=20;      // number of satelite circles
const double R=10.0;    // satelite circles radius
const double E= 0.7;    // ellipse distortion ry=rx*E
struct _circle { double x,y,r; _circle() { x=0; y=0; r=0.0; } } circle[N];
struct _ellipse { double x,y,rx,ry; _ellipse() { x=0; y=0; rx=0.0; ry=0.0; } } ellipse;

int i,j,k;
double l,a,da,m,dm,x,y,q,r0;
l=double(N)*R;                          // circle cener lines polygon length
ellipse.x =0.0;                         // set ellipse parameters
ellipse.y =0.0;
r0=divide(l,M_PI*sqrt(0.5*(1.0+(E*E))))-R;// aprox radius to match ellipse length for start
l=R+R; l*=l;
m=1.0; dm=1.0; x=0.0;
for (k=0;k<5;k++)                       // aproximate ellipse size to the right size
    {
    dm=fabs(0.1*dm);                    // each k-iteration layer is 10x times more accurate
    if (x>l) dm=-dm;
    for (;;)
        {
        ellipse.rx=r0  *m;
        ellipse.ry=r0*E*m;
        for (a=0.0,i=0;i<N;i++)         // set circle parameters
            {
            q=(2.0*a)-atanxy(cos(a),sin(a)*E);
            circle[i].x=ellipse.x+(ellipse.rx*cos(a))+(R*cos(q));
            circle[i].y=ellipse.y+(ellipse.ry*sin(a))+(R*sin(q));
            circle[i].r=R;    
            da=divide(360*deg,N); a+=da;
            for (j=0;j<5;j++)           // aproximate next position to match 2R distance from current position
                {
                da=fabs(0.1*da);        // each j-iteration layer is 10x times more accurate
                q=(2.0*a)-atanxy(cos(a),sin(a)*E);
                x=ellipse.x+(ellipse.rx*cos(a))+(R*cos(q))-circle[i].x; x*=x;
                y=ellipse.y+(ellipse.ry*sin(a))+(R*sin(q))-circle[i].y; y*=y; x+=y;
                if (x>l) for (;;)       // if too far dec angle
                    {
                    a-=da;
                    q=(2.0*a)-atanxy(cos(a),sin(a)*E);
                    x=ellipse.x+(ellipse.rx*cos(a))+(R*cos(q))-circle[i].x; x*=x;
                    y=ellipse.y+(ellipse.ry*sin(a))+(R*sin(q))-circle[i].y; y*=y; x+=y;
                    if (x<=l) break;
                    }
                else if (x<l) for (;;)  // if too short inc angle
                    {
                    a+=da;
                    q=(2.0*a)-atanxy(cos(a),sin(a)*E);
                    x=ellipse.x+(ellipse.rx*cos(a))+(R*cos(q))-circle[i].x; x*=x;
                    y=ellipse.y+(ellipse.ry*sin(a))+(R*sin(q))-circle[i].y; y*=y; x+=y;
                    if (x>=l) break;
                    }
                else break;
                }
            }
        // check if last circle is joined as it should be
        x=circle[N-1].x-circle[0].x; x*=x;
        y=circle[N-1].y-circle[0].y; y*=y; x+=y;
        if (dm>0.0) { if (x>=l) break; }
        else        { if (x<=l) break; }
        m+=dm;
        }
    }

这段代码可能有些混乱，以下是一些信息：

1. 首先尝试将椭圆的rx、ry轴设置得尽可能接近

   椭圆的长度应该约为N*R*2，即圆心之间线段的多边形长度。

2. 尝试组合圆，使它们彼此相连并与椭圆相接触

   我使用椭圆角度的迭代来实现。问题在于圆在其位置角度上不与椭圆接触，因此有q变量...以在椭圆法线周围进行补偿。请查看图像中黄色金色的线。

3. 放置圆后检查最后一个圆是否与第一个圆接触

   如果没有，则插值椭圆的大小，实际上通过m变量向上或向下缩放rx,ry。

4. 您可以调整精度

   通过更改j,k for和/或更改dm,da缩放因子。

5. 输入参数E应至少为0.5，最大为1.0

   否则，圆的位置可能会错位，因为在非常偏心的椭圆上不可能放置圆（如果N太低）。理想设置是0.7<=E<=1.0，越靠近1算法越安全。

6. atanxy(dx,dy)与`atan(dy/dx)相同

   但它通过对dx,dy的符号分析处理所有4个象限，就像atan2(dy,dx)一样。

希望能对您有所帮助。

- Spektre

---

Spektre，谢谢你理解我想要和我知道的内容。然而，在圆形情况下，我知道如何获取每个卫星的x和y。这种情况很容易，因为每个卫星的中心在内部圆周围形成一个正多边形。但是在椭圆情况下，由于每个角度都不同，所以无法使用cos和sin来找到x和y。 - user2360915

太好了。非常感谢。代码有点凌乱，但我知道了大概意思。谢谢。 - user2360915

如果您需要使用E<0.7，则应更改q计算（首先尝试放置2个圆，然后使用前一个和下一个圆之间的中心线之间的中间角度），这对于较小的E更精确，因为较小的E具有更高的精度损失。 - Spektre