2013-04-17 78 views
1

我有一個粒子系統,我使用boost::geometry將我的橢圓粒子近似爲多邊形,然後使用庫的交集函數查找重疊區域。我正在計算一個「內部」和「外部」橢圓(多邊形)區域,爲每個粒子 - 粒子相互作用分配一個「潛在」。具有Boost :: geometry幾何嚴重性能下降的多邊形交叉點

我的潛在功能是這樣的:

double Potential(Cell* current, Cell* next) 
{ 
    double areaRep, areaAtt; 
    double distance = Distance(current,next); 
    double A1 = current->getLength(); 
    double B1 = A1/2.0; 
    double theta1 = current->getTheta(); //*180.0/M_PI 
    double x1 = current->getCurrX(); 
    double y1 = current->getCurrY(); 
    double A2 = next->getLength(); 
    double B2 = A2/2.0; 
    double theta2 = next->getTheta(); 
    double x2 = next->getCurrX(); 
    double y2 = next->getCurrY(); 
    polygon_2d poly1, poly2, poly3, poly4; 
    double lamda1, lamda2; 
    lamda1 = 0.0005; lamda2 = 0.00001; 
    if(distance < 2.0*1.5*A1) { 
     ellipse2poly(theta1, A1, B1, x1, y1, &poly1); 
     ellipse2poly(theta2, A2, B2, x2, y2, &poly2); 
     areaRep = getOverlapingAreaPoly(poly1,poly2); 
     ellipse2poly(theta1, 1.5*A1, 1.5*B1, x1, y1, &poly3); 
     ellipse2poly(theta2, 1.5*A2, 1.5*B2, x2, y2, &poly4); 
     areaAtt = getOverlapingAreaPoly(poly3, poly4); 
     return (lamda1*areaRep - lamda2*areaAtt); 
    } 
    else 
     return 0.0; 
} 

的 「多邊形化」 功能是:

int ellipse2poly(double theta, double A1, double B1, double H1, double K1, polygon_2d *po) 
{ 
    using namespace boost::geometry; 
    polygon_2d poly; 
    const int n = 20; 
    double angle = theta; // cell orientation 
    double a = A1; // Long semi-axis length 
    double b = B1; // short semi-axis length 
    double xc = H1; // current X position 
    double yc = K1; // current Y position 

    if(!n) 
    { 
     std::cout << "error ellipse(): n should be >0\n" <<std::endl; 
     return 0; 
    } 
    double t = 0; 
    int i = 0; 
    double coor[2*n+1][2]; 

    double x, y; 
    double step = M_PI/(double)n; 
    double sinphi = sin(angle); 
    double cosphi = cos(angle); 
    for(i=0; i<2*n+1; i++) 
    { 
     x = xc + a*cos(t)*cosphi - b*sin(t)*sinphi; 
     y = yc + a*cos(t)*sinphi + b*sin(t)*cosphi; 
     coor[i][0] = x; 
     coor[i][1] = y; 
     t += step; 
    } 
    assign_points(poly, coor); 
    correct(poly); 
    *po = poly; 
    return 1; 
} 

並且返回的區域是:

double getOverlapingAreaPoly(polygon_2d poly, polygon_2d poly2) 
{ 
    point_2d cent; //centre of overlaping area 
    double overAreaPoly = 0.0; 
    typedef std::vector<polygon_2d > polygon_list; 
    polygon_list v; 

    intersection(poly,poly2,v); 
    for (polygon_list::const_iterator it = v.begin(); it != v.end(); ++it) 
    { 
     centroid(*it, cent); 
     overAreaPoly = area(*it); 
    } 

    return overAreaPoly; 
} 

的函數被調用每一個細胞(粒子),只要它不是同一個細胞。以前,使用另一種方法,對於100個粒子的一次迭代,我的算法的一次迭代大約需要43 ms。現在大約需要1分鐘(!!!),所以我想我已經做了一些可怕的錯誤!

我已經測試過這隻在MSVC2012下win7 64位。我將用Qt 4.7.4報告Linux Mint。

編輯: 我已經在Linux Mint上用Qt 4.7.4進行了測試,它運行非常合理;也許90-100毫秒每迭代這很好。我不知道win7中有什麼問題...

+0

是你的潛力()函數是固定的還是你願意修改它?你能舉一個例子說明你打算如何使用潛力? – Matthias

+1

其實我修好了。我剛剛在Visual Studio中啓動了一個新項目,複製了源文件和頭文件,重新編譯,並且它再次合理地工作!我不知道爲什麼它實際上是... 爲了回答你的問題,我使用粒子之間距離的歸一化矢量來確定相互作用的方向,然後我將它與勢能相乘。 – Dima1982

+0

謝謝你的澄清。很高興聽到它現在工作正常:-) – Matthias

回答

1

我已經修復了它。我在Visual Studio中開始了一個新項目,並複製了所有的源文件和頭文件,重新編譯,現在一切運行順利。我想從根本上改變代碼和加/減東西必須有一些影響...

相關問題