00001
00002
00003
00004 #include "tensor_field_nav_core/QuadMesh.h"
00005
00006 QuadMesh::QuadMesh()
00007 {
00008 quad_verts = NULL;
00009 quadcells = NULL;
00010 elist = NULL;
00011 edgelist = NULL;
00012
00013 nverts = nfaces = nedges = 0;
00014 }
00015
00016
00017 QuadMesh::QuadMesh(int xdim, int ydim, double xstart, double xend,
00018 double ystart, double yend)
00019 {
00020 quad_verts = NULL;
00021 quadcells = NULL;
00022 elist = NULL;
00023 edgelist = NULL;
00024
00025 nverts = nfaces = nedges = 0;
00026
00027 gen_regquad_mesh(xdim, ydim, xstart, xend, ystart, yend);
00028
00029
00030 construct_edges();
00031 orient_edges_cells();
00032
00033 XDIM = xdim;
00034 YDIM = ydim;
00035 this->xstart = xstart;
00036 this->xend = xend;
00037 this->ystart = ystart;
00038 this->yend = yend;
00039 }
00040 QuadMesh::~QuadMesh()
00041 {
00042 finalize_quad_verts();
00043 finalize_quad_cells();
00044 }
00045
00046
00047 void QuadMesh::finalize_quad_verts()
00048 {
00049 int i;
00050 if(quad_verts != NULL)
00051 {
00052 for(i=0; i<nverts; i++)
00053 {
00054 if(quad_verts[i]!=NULL)
00055 {
00056 if(quad_verts[i]->edges != NULL)
00057 free(quad_verts[i]->edges);
00058 free(quad_verts[i]);
00059 }
00060 }
00061 free(quad_verts);
00062 }
00063 }
00064
00065
00066 bool QuadMesh::gen_regquad_mesh(int xdim, int ydim, double xstart, double xend,
00067 double ystart, double yend)
00068 {
00069
00070 if(gen_regquad_vertices(xdim, ydim, xstart, xend, ystart, yend))
00071 {
00072
00073 gen_regquad_faces(xdim, ydim);
00074
00075 return true;
00076 }
00077 else
00078 return false;
00079 }
00080
00081
00082
00083 bool QuadMesh::gen_regquad_vertices(int xdim, int ydim, double xstart, double xend,
00084 double ystart, double yend)
00085 {
00086 if(xstart > xend || yend < ystart) return false;
00087 int i, j;
00088 double d_x, d_y, x_coord, y_coord;
00089 xinterval = d_x = (xend-xstart)/(xdim-1);
00090 yinterval = d_y = (yend-ystart)/(ydim-1);
00091 int quad_vert_id;
00092
00093
00094 if(quad_verts != NULL)
00095 finalize_quad_verts();
00096 nverts = xdim*ydim;
00097 quad_verts = (QuadVertex**)malloc(sizeof(QuadVertex*)*nverts);
00098 for(i=0; i<nverts; i++)
00099 quad_verts[i] = (QuadVertex*)malloc(sizeof(QuadVertex));
00100
00101 for(j=0; j<ydim; j++)
00102 {
00103 y_coord = ystart+j*d_y;
00104 for(i=0; i<xdim; i++)
00105 {
00106 quad_vert_id = j*xdim+i;
00107 x_coord = xstart+i*d_x;
00108 quad_verts[quad_vert_id]->index = quad_vert_id;
00109 quad_verts[quad_vert_id]->x = x_coord;
00110 quad_verts[quad_vert_id]->y = y_coord;
00111 }
00112 }
00113
00114 init_vertices();
00115
00116 return true;
00117 }
00118
00119
00120
00121 void QuadMesh::init_vertices()
00122 {
00123 int i;
00124 QuadVertex *qv;
00125 for(i=0; i<nverts; i++)
00126 {
00127 qv = quad_verts[i];
00128 qv->edges = NULL;
00129 qv->Num_edge = 0;
00130 qv->InRegion = false;
00131 qv->distance = 1.e50;
00132 qv->OnBoundary = false;
00133 qv->RegionListID = -1;
00134
00135 qv->visited = false;
00136 qv->Jacobian.set(0.);
00137 qv->major.set(0.);
00138 qv->minor.set(0.);
00139 qv->ncells = 0;
00140 qv->cells = NULL;
00141
00142 qv->inland = true;
00143 qv->inveg = false;
00144
00145 qv->which_region=0;
00146 qv->inbrushregion=false;
00147
00148 qv->phi=0.;
00149
00150 qv->density=1.;
00151 }
00152 }
00153
00154
00155
00156
00157 void QuadMesh::finalize_quad_cells()
00158 {
00159 int i;
00160 if(quadcells != NULL)
00161 {
00162 for(i=0; i<nfaces; i++)
00163 {
00164 if(quadcells[i] != NULL)
00165 {
00166 free(quadcells[i]->verts);
00167 free(quadcells[i]);
00168 }
00169 }
00170 free(quadcells);
00171 quadcells=NULL;
00172 }
00173 }
00174
00175
00176
00177
00178
00179 void QuadMesh::gen_regquad_faces(int xdim, int ydim)
00180 {
00181 int i, j;
00182 QuadCell *f;
00183
00184
00185 if(quadcells != NULL)
00186 {
00187 finalize_quad_cells();
00188 }
00189
00190 nfaces = (xdim-1)*(ydim-1);
00191 quadcells = (QuadCell**)malloc(sizeof(QuadCell*)*nfaces);
00192 for(i=0; i<nfaces; i++)
00193 {
00194 quadcells[i] = (QuadCell*)malloc(sizeof(QuadCell));
00195 quadcells[i]->verts = (int*)malloc(sizeof(int)*4);
00196 }
00197
00198 int cell_index = 0;
00199
00200 for(j=0; j<ydim-1; j++)
00201 {
00202 for(i=0; i<xdim-1; i++)
00203 {
00204 cell_index = j*(xdim-1)+i;
00205 quadcells[cell_index]->verts[0] = j*xdim+i;
00206 quadcells[cell_index]->verts[1] = j*xdim+(i+1);
00207 quadcells[cell_index]->verts[2] = (j+1)*xdim+(i+1);
00208 quadcells[cell_index]->verts[3] = (j+1)*xdim+i;
00209 quadcells[cell_index]->nverts=4;
00210 quadcells[cell_index]->index = cell_index;
00211
00212
00213 quadcells[cell_index]->xstart = i;
00214 quadcells[cell_index]->xend = i+1;
00215 quadcells[cell_index]->ystart = j;
00216 quadcells[cell_index]->yend = j+1;
00217 quadcells[cell_index]->x_start_coord = quad_verts[quadcells[cell_index]->verts[0]]->x;
00218 quadcells[cell_index]->y_start_coord = quad_verts[quadcells[cell_index]->verts[0]]->y;
00219
00220 quadcells[cell_index]->degpt_index = -1;
00221
00222 quadcells[cell_index]->OnBoundary=false;
00223 quadcells[cell_index]->visited=false;
00224
00225
00226
00227
00228
00229
00230 quadcells[cell_index]->maj_samplepts=NULL;
00231 quadcells[cell_index]->maj_nsamplepts = 0;
00232 quadcells[cell_index]->MAJMaxSampNum = 0;
00233 quadcells[cell_index]->min_samplepts=NULL;
00234 quadcells[cell_index]->min_nsamplepts = 0;
00235 quadcells[cell_index]->MINMaxSampNum = 0;
00236
00237
00238 quadcells[cell_index]->majorlines = NULL;
00239 quadcells[cell_index]->hasmajor = false;
00240 quadcells[cell_index]->minorlines = NULL;
00241 quadcells[cell_index]->hasminor = false;
00242
00243 quadcells[cell_index]->sketchlines = NULL;
00244
00245 quadcells[cell_index]->intersectlist = NULL;
00246 quadcells[cell_index]->nintersects = 0;
00247
00248 quadcells[cell_index]->streetgraphedgelist = NULL;
00249 quadcells[cell_index]->nstreetgraphedges = 0;
00250
00251 quadcells[cell_index]->which_region = 0;
00252 quadcells[cell_index]->in_region=false;
00253
00254 quadcells[cell_index]->mapbounds=NULL;
00255 quadcells[cell_index]->nmapbounds=0;
00256
00257
00258
00259
00260
00261 QuadVertex *v;
00262 v = quad_verts[quadcells[cell_index]->verts[0]];
00263 v->cells= extend_celllist_ver(v->cells, v->ncells);
00264 v->cells[v->ncells] = quadcells[cell_index];
00265 v->ncells ++;
00266
00267 v = quad_verts[quadcells[cell_index]->verts[1]];
00268 v->cells= extend_celllist_ver(v->cells, v->ncells);
00269 v->cells[v->ncells] = quadcells[cell_index];
00270 v->ncells ++;
00271
00272 v = quad_verts[quadcells[cell_index]->verts[2]];
00273 v->cells= extend_celllist_ver(v->cells, v->ncells);
00274 v->cells[v->ncells] = quadcells[cell_index];
00275 v->ncells ++;
00276
00277 v = quad_verts[quadcells[cell_index]->verts[3]];
00278 v->cells= extend_celllist_ver(v->cells, v->ncells);
00279 v->cells[v->ncells] = quadcells[cell_index];
00280 v->ncells ++;
00281 }
00282 }
00283 }
00284
00285 void QuadMesh::construct_edges()
00286 {
00288 QuadEdge *Cur_elink;
00289
00290 elist = new QuadEdge();
00291 elist->index = -1;
00292 elist->next = NULL;
00293 Cur_elink = elist;
00294
00295 int edge_id = 0;
00296 nedges = 0;
00297
00299 QuadVertex **vert = quad_verts;
00300
00302 int i, j, m, n;
00303 int Cur_vert, Next_vert;
00304
00305 for( i = 0; i < nfaces; i++)
00306 {
00307 for( j = 0; j < quadcells[i]->nverts; j++)
00308 {
00309
00310 if( j == quadcells[i]->nverts - 1){
00311 Cur_vert = quadcells[i]->verts[j];
00312 Next_vert = quadcells[i]->verts[0];
00313 }
00314 else{
00315 Cur_vert = quadcells[i]->verts[j];
00316 Next_vert = quadcells[i]->verts[j+1];
00317 }
00318
00319
00320 if(vert[Cur_vert]->Num_edge == 0 || vert[Next_vert]->Num_edge == 0)
00321
00322 {
00324 QuadEdge *new_edge = new QuadEdge;
00325 new_edge->index = edge_id;
00326 edge_id ++;
00327
00329 new_edge->tris[0] = -1;
00330 new_edge->tris[1] = -1;
00331
00332 new_edge->tris[0] = i;
00333 new_edge->verts[0] = Cur_vert;
00334 new_edge->verts[1] = Next_vert;
00335 new_edge->visited = false;
00336 new_edge->next = NULL;
00337 Cur_elink->next = new_edge;
00338 Cur_elink = new_edge;
00339
00340
00341
00342
00343
00344 vert[Cur_vert]->edges =
00345 Extend_Elist(vert[Cur_vert]->edges, vert[Cur_vert]->Num_edge);
00346 vert[Next_vert]->edges =
00347 Extend_Elist(vert[Next_vert]->edges, vert[Next_vert]->Num_edge);
00348
00349 vert[Cur_vert]->edges[vert[Cur_vert]->Num_edge] = new_edge;
00350 vert[Next_vert]->edges[vert[Next_vert]->Num_edge] = new_edge;
00351
00352 vert[Cur_vert]->Num_edge++;
00353 vert[Next_vert]->Num_edge++;
00354
00356 quadcells[i]->edges[j] = new_edge;
00357
00359 nedges++;
00360 }
00361 else{
00362 for( m = 0; m < vert[Cur_vert]->Num_edge; m++)
00363 for( n = 0; n < vert[Next_vert]->Num_edge; n++)
00364 {
00365 if( vert[Cur_vert]->edges[m]->index
00366 == vert[Next_vert]->edges[n]->index)
00367
00368 {
00369
00370 vert[Cur_vert]->edges[m]->tris[1] = i;
00371
00372 quadcells[i]->edges[j] = vert[Cur_vert]->edges[m];
00373
00374 goto LL;
00375
00376 }
00377 }
00378
00379 LL: if( m > vert[Cur_vert]->Num_edge - 1 )
00380
00381 {
00383 QuadEdge *new_edge = new QuadEdge;
00384 new_edge->index = edge_id;
00385 edge_id ++;
00387 new_edge->tris[0] = -1;
00388 new_edge->tris[1] = -1;
00389
00390 new_edge->tris[0] = i;
00391 new_edge->verts[0] = Cur_vert;
00392 new_edge->verts[1] = Next_vert;
00393 new_edge->visited = 0;
00394 new_edge->next = NULL;
00395 Cur_elink->next = new_edge;
00396 Cur_elink = new_edge;
00397
00398
00399
00400
00402
00403 vert[Cur_vert]->edges =
00404 Extend_Elist(vert[Cur_vert]->edges, vert[Cur_vert]->Num_edge);
00405 vert[Next_vert]->edges =
00406 Extend_Elist(vert[Next_vert]->edges, vert[Next_vert]->Num_edge);
00407
00408 vert[Cur_vert]->edges[vert[Cur_vert]->Num_edge] = new_edge;
00409 vert[Next_vert]->edges[vert[Next_vert]->Num_edge] = new_edge;
00410
00411 vert[Cur_vert]->Num_edge++;
00412 vert[Next_vert]->Num_edge++;
00413
00415
00416 quadcells[i]->edges[j] = new_edge;
00417
00419 nedges++;
00420 }
00421 }
00422
00423 }
00424 }
00425
00426
00427
00428
00429
00430 edgelist=(QuadEdge**)malloc(sizeof(QuadEdge*)*nedges);
00431 QuadEdge *e = elist->next;
00432 for(i=0; i<nedges; i++)
00433 {
00434 edgelist[i]= e;
00435 e=e->next;
00436 }
00437 }
00438
00439
00440 QuadEdge **QuadMesh::Extend_Elist(QuadEdge **edge_link, int Num_edges)
00441 {
00442 QuadEdge **temp = edge_link;
00443 QuadEdge **new_edge_link = (QuadEdge **) malloc (sizeof (QuadEdge*)*(Num_edges+1));
00444 if( Num_edges > 0)
00445 {
00446 for(int i = 0; i < Num_edges; i++)
00447 new_edge_link[i] = temp[i];
00448 free (temp);
00449 }
00450
00451 return new_edge_link;
00452 }
00453
00454
00455 QuadCell **QuadMesh::extend_celllist_ver(QuadCell **cells, int ncells)
00456 {
00457 QuadCell **temp = cells;
00458 QuadCell **new_cells = (QuadCell **) malloc (sizeof (QuadCell*)*(ncells+1));
00459 if( ncells > 0)
00460 {
00461 for(int i = 0; i < ncells; i++)
00462 new_cells[i] = temp[i];
00463 free (temp);
00464 }
00465 cells = new_cells;
00466 return new_cells;
00467 }
00468
00469
00470
00471
00472
00473
00474
00475
00476 void QuadMesh::orient_edges_cells()
00477 {
00478 int i, j, k;
00479 QuadCell *face;
00480 QuadEdge *e[4];
00481
00482 for(i=0; i<nfaces; i++)
00483 {
00484 face = quadcells[i];
00485 for(j=0; j<face->nverts; j++)
00486 {
00487 if((face->edges[j]->verts[0]==face->verts[0]&&face->edges[j]->verts[1]==face->verts[1])
00488 ||(face->edges[j]->verts[0]==face->verts[1]&&face->edges[j]->verts[1]==face->verts[0]))
00489 e[0] = face->edges[j];
00490 else if((face->edges[j]->verts[0]==face->verts[1]&&face->edges[j]->verts[1]==face->verts[2])
00491 ||(face->edges[j]->verts[0]==face->verts[2]&&face->edges[j]->verts[1]==face->verts[1]))
00492 e[1] = face->edges[j];
00493 else if((face->edges[j]->verts[0]==face->verts[2]&&face->edges[j]->verts[1]==face->verts[3])
00494 ||(face->edges[j]->verts[0]==face->verts[3]&&face->edges[j]->verts[1]==face->verts[2]))
00495 e[2] = face->edges[j];
00496 else
00497 e[3] = face->edges[j];
00498 }
00499
00500 for(j=0; j<face->nverts; j++)
00501 face->edges[j]=e[j];
00502 }
00503 }