CGAL users discussion list

[parisa rahmani <>]

Hi,

I want to simulate a system of particles for N steps. In every step, periodic triangulation of these points is constructed. In every step coordinates of particles change, So, I need to remove (or delete) triangulation of the previous step and make a new one.

In step 0 every thing is good, but in the next step, when triangulation is going to be remade from modified points, SEGMENTATION FAULT appears.

Then according to suggestion of someone from Stackoverflow.com, I saved the new points (points updated in step 0) in a file with double precision and tried to construct a triangulation from these points,  in this case too, the SEG-FAULT appears.

using namespace std;

//******************CGAL_prerequisite************************/
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef CGAL::Periodic_2_triangulation_traits_2<Kernel> Gt;
typedef CGAL::Triangulation_vertex_base_with_info_2<unsigned int, Gt> Vb;
typedef CGAL::Periodic_2_triangulation_face_base_2<Gt> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> Tds;
typedef CGAL::Periodic_2_Delaunay_triangulation_2<Gt, Tds> Triangulation;

typedef Triangulation::Point Point;
typedef Triangulation::Edge_iterator Edge_iterator;
typedef Triangulation::Vertex_handle Vertex_handle;
typedef Triangulation::Vertex_circulator Vertex_circulator;
typedef Triangulation::Iso_rectangle Iso_rectangle;
typedef Triangulation::Vertex_iterator Vi;

typedef std::vector<std::pair<Point, unsigned> > Vector;

/*********************Simulation_parameters*****************/
const int N = 16;            
double L = 4;        //linear size of the square 
double S = 0.5;         //0.5
double dt = 1.0;        
int N_step = 100;

ofstream output("traj_2D.xyz");

int main(){

Iso_rectangle domain(0,0,L,L);

Vector points;

vector<vector<double> > dir;//an array for 2D velocities of each particle
dir.resize(N);
for(int i=0; i<N; i++){dir[i].resize(2);}

vector<vector<double> > pos;//an array for 2D initial positions
pos.resize(N);
for(int i=0; i<N; i++){pos[i].resize(2);}


/*****initialize_positions*****/
//some part of the code
/*****initialize_velocities*****/
//some part of the code   


  for (int i = 0; i < N; i++)
       points.push_back(make_pair(Point(pos[i][0], pos[i][1]), i));  

/********************************************************************
*********************************************************************/

for(int istep = 0; istep < N_step; istep++)
     {
      Triangulation T(points.begin(), points.end(), domain); 


      //update positions and apply periodic boundary conditons
        for (int i = 0; i < N; i++) {

             points[i].first = Point(points[i].first.x()+S*dir[i][0]*dt , points[i].first.y()+S*dir[i][1]*dt);


            // use periodic boundary conditions
            if (points[i].first.x() < 0)
                points[i].first=Point(points[i].first.x()+L, points[i].first.y());
            if (points[i].first.y() < 0)
                points[i].first=Point(points[i].first.x(), points[i].first.y()+L);
            if (points[i].first.x() > L)
                points[i].first=Point(points[i].first.x()-L, points[i].first.y());
            if (points[i].first.y() > L)
                points[i].first=Point(points[i].first.x(), points[i].first.y()+L);

                                 }
          //write the updated positions to a file

  output.precision(17);
  for(int s = 0; s < N; s++){
             output<<points[s].first.x()<<"  "<<points[s].first.y()<<"  "<<"0.0"<<endl;}
  T.clear();//I don't know, whether this is necessary or not.
  ifstream in("traj_2D.xyz"); 
  istream_iterator<Point> begin(in);
  istream_iterator<Point> end;
  Triangulation TRI;
 TRI.insert(begin, end);}    
return 0;}

Any suggestion would be appreciated.

Parisa